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Abstract 


Dale, Virginia H.; Hemstrom, Miles. CLIMACS: A computer model of forest stand 
development for western Oregon and Washington. Res. Pap. PNW-327. 
Portland, OR: U.S. Department of Agriculture, Forest Service, Pacific Northwest 
Forest and Range Experiment Station; 1984. 60 p. 


A simulation model for the development of timber stands in the Pacific Northwest 
is described. The model grows individual trees of 21 species in a 0.20-hectare 
(0.08-acre) forest gap. The model provides a means of assimilating existing infor- 
mation, indicates where knowledge is deficient, suggests where the forest system 
is most sensitive, and provides a first testing ground for hypotheses. Model verifi- 
cation simulations are included for up to 500 years on various sites. Fire, wind, or 
clearcutting can occur at intervals and intensities specified by users. The model 
was developed by modifying an existing forest succession simulator of eastern 
deciduous forests. Birth, growth, and death of individual trees are functions of 
existing light and temperature conditions, competition and species characteristics. 
Modifications of the existing simulator include tree height growth being related 
to temperature and moisture conditions, the foliage biomass to diameter relation- 
ship being more realistic, and five mortality classes and shade tolerance classes 
being defined. 


Keywords: Succession, models, simulation, community dynamics (plant). 
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Introduction 


CLIMACS (Computer Linked Integrative Model for Assessing Community Struc- 
ture) is a simulation model of forest succession for western Oregon and Washing- 
ton. The model tracks characteristics of individual trees growing in a forest open- 
ing of 0.20 hectare (0.08 acre). The FORTRAN IV code of CLIMACS is based on 
FORET (Shugart and West 1977), a succession model for eastern Tennessee, 
which in turn was modified from JABOWA (Botkin and others 1972), a simulator 
for forests of the northeastern United States. CLIMACS retains the stochastic 
features of the succession models for the eastern deciduous forests and has a 
more elaborate diameter increment equation and treatment of mortality. A mois- 
ture stress index is incorporated into the function that calculates diameter incre- 
ment, and height growth is dependent upon site quality for some species. Mortality 
is related to the size of the tree and the successional status of the species under con- 
sideration. Three types of disturbance can occur: fire, wind, or clearcuts. 


The model paradigm represents our view of succession as the result of individual 
tree processes and episodic disturbances. The term “climax” is tenuous because 
pioneer species are so long-lived that climax forests rarely develop in the Pacific 
Northwest even though major disturbances are infrequent. Although the simulator 
should be realistic enough to approximate stand development on particular sites, 
its main value is in summarizing and testing current hypotheses about tree growth 
and the factors affecting succession. The model is an exercise that should be 
useful to managers to examine hypotheses of the long-term effects of human 
activities on stand development. 


A review of forest succession models is presented by Shugart and West (1980). 
According to their terminology, CLIMACS is a gap model that simulates charac- 
teristics of individual trees in a small portion of the forest—either an opening in 
the canopy or a sample plot. Such models are particularly useful for evaluation of 
long-term and large-scale changes in the environment and the effects of those 
changes on forest succession. 


The purpose of this paper is to clarify algorithms and subroutines of the model. 
The basic features of the model are described in Botkin and others (1972) and 
Shugart and West (1977). Following the standards for model documentation 

from Swartzman (1979), the assumptions of the model are set forth, the value and 
source of each parameter is given, the source of the equations and potential 
problems of the curve forms are indicated, and the logic of the FORTRAN pro- 
gram is presented. This paper describes the structure of the main program and 
discusses the four major subroutines: BIRTH, GROW, KILL, and DISTRB. Results 
for several forest sites are presented. 


The model code is stored on computer tape at the Forest Science Data Bank, 
Department of Forest Science, Oregon State University, Corvallis, Oregon. 


Structure and 
Assumptions 


Table 1—State variables’ 


Name Definition Units Subroutines 

A(KYR,i) Aboveground biomass for species t/ha MAIN, OUTPUT, PLTREE 

A(KYR, NSPEC+2) Total number of trees per plot 

A(KYR, NSPEC+3) Leaf area index for plot m2/m2 PLTREE, OUTPUT, MAIN 

A(KYR, NSPEC+1) Stand biomass t/ha 

DBH(j) Diameter at breast height cm SPROUT, KILL, GROW, OUT- 
= of tree j PUT, BIRTH 

IAGE(k) Age of tree k yrs GROW, KILL, SPROUT 

NTREES(i) Number of trees of species i SPROUT, GROW 

BIOM Stand biomass = A(KYR,i) t/ha PLTREE 

BAREA Basal area of stand — m¢“/ha PLTREE 

A(KYR1, 26) Total foliage biomass t/ha OUTPUT, PLTREE 

FB(i) Foliage biomass for species t/ha OUTPUT, PLTREE 

D(i,k) Number of stems of species i PLTREE 


in diameter class (k*10) to 


(k*10 + 10) 


INumbers that represent the current state or condition of the system. 


The model updates and records on an annual basis characteristics of all individ- 
ual trees greater than 10 cm diameter at breast height (d.b.h.) on a 0.20-hectare 
(0.08-acre) plot (table 1). In the model, foliage is distributed vertically but not 
horizontally (Cartesian coordinates of each tree are not recorded). The driving 
variables (numbers that are inputs to the model and are not affected by the com- 
ponents of the system) for the model are plant moisture stress (PMS) (-bars) and 
temperature growth index (TGII) (days). PMS is the predawn negative xylem pres- 
sure measured near the end of the growing season (Waring and Cleary 1967). 
TGIlI is from a temperature summation formula that weighs temperatures by their 
effects on the production of seedlings of Pseudotsuga menziesii (Mirb.) Franco 
(Cleary and Waring 1969). In the model, the yearly temperature growth index 
(TGI) varies from year to year in a normal distribution around the driving variable 
TGIIl. The model has six major assumptions: 


1. PMS and TGII reflect the main physical driving variables of forest succession. 

2. Competition for light is the primary biological factor affecting forest succession. 
3. Regeneration of species that can grow under existing light and soil conditions is 
stochastic within the bounds of the geographic distribution of the species. 

4. Diameter increment can be modeled as a multiplicative function of tree size, photo- 
synthetic rate as reflected in foliage biomass, nutrient competition, light availability, 
shade tolerance, and moisture and temperature stress. 

5. Mortality can occur at two levels—an individual tree and a portion of the stand. 
6. Individual tree mortality is conditioned by growth rate. 


A listing of the entire program is given in Appendix 1. 


The model can be applied to four geographic regions west of the Cascade Range: 
Santiam Pass south to the California border, Santiam Pass north to Mount 
Rainier, north of Mount Rainier to the Canadian border, and the Olympic Penin- 
sula (fig. 1). Restriction of a model run to a geographic region limits the species 
in the succession model to those found naturally in the region. 


Main Program 


Mount }jRainier 


Helens 
nivel 
Columbia 
Santiam Pass 


OREGON 


Figure 1.—The model is applicable to four geographic regions 
in Oregon and Washington. 


The main program manages the entire model by initializing the year and plot rep- 
licate counters, setting the data storage array to zero, and calling the major sub- 
routines as outlined in figure 2. A successional sequence can occur on replicate 
plots and for as many years as specified by the user. The model can start from a 
bare plot (see subroutine BIRTH), or the user can read in the species and diam- 
eter of each tree on the plot (see subroutine INPUT). The names and parameters 
of the species, driving variables, and length of the model run are read in by a call 
from the main program to subroutine INPUT. For each year of the simulation the 
appropriate species for the site are selected (Subroutine SELECT), trees can die 
(subroutine KILL), regeneration occurs (Subroutine BIRTH), and the diameter of 


INITIALIZE model parameters 


Call SELECT: choose species 
appropriate for geographic region 
Call KILL: mortality occurs | 
Call BIRTH 


Call GROW: increment 


d.b.h. of trees Cycle 
ry for 
specified 
Call OUTPUT: records total number 
biomass, number of trees and 
leaf area index of 
years 


Catastrophic 
disturbance 


Call DISTRB 


Figure 2.—General scheme of model. 


Table 2—Subroutines 


Subroutine 


BIRTH(KYR) 


ERR 
GGNORD 
(NSEED1, NSEED2, Z) 
GROW 
INIT 


INPUT 


KILL 


OUTPUT 
PLOTIN(IPLOT ) 


PLTREE 
RANDOM( NSEED) 
RANSEE 


SEERECT 


SPROUT 


Discussion 


Introduces new trees with the diameter randomly 
distributed between 10 and 15 cm at breast height 
of the tree. KYR = current year. 


Indicates error when there are more than 700 
trees per plot. 


Provides a normal distribution, Z, based on 
the random numbers NSEED1 and NSEED2. 


Growth of d.b.h. and leaf area. 


NOGRO(I) 
KSPRT(1) 


OO! 
1, NSPEC. 


OF 
eal 


Reads and writes species specific names and 
parameters. 


Eliminates trees as a function of diameter and 
growth rate. 


Updates biomass and the leaf area index. 


Variables to start simulation in a bare plot 
(IPLOT = plot number). 


Writes out information on species and years. 
Random number generator (NSEED = seed). 


Reads random number generator seed from input 
file. 


Checks water stress and temperature requirements 
for species in the four geographic regions. 


Stump sprouting for appropriate species. 


existing trees increases (Subroutine GROW). A discussion of each subroutine is 
given in table 2. For specified years, stand and species characteristics are printed 
out in tabular form (see Appendix 2). 


Definitions of the constants, indices, and parameters are given in tables 3, 4, 5, 
and 6. The parameter values are appropriate for species in the Pacific Northwest 
and their sources are indicated in tables 7, 8 and 9. 


Table 3—Constants 


Name 


soilq 
size 


PI2 
phi 


tinc 


Value 


100000.0 


Definition Units 


Maximum biomass for kg/ha 
forests in the area 
Mean d.b.h. of new trees cm 
27 

Proportion of annual 

light attenuation 

constant 
Minimum growth of d.b.h. cm 
per year to trigger 


mortality 
Table 4—Indexes 
Index Range Meaning 
ab 1 - NSPEC Species of tree 
i 1 - NTOT Individual tree 
KYR 1 - NYEAR Model year 
IPLOT 1 - total number Plot number 


of plots 


Subroutines 


GROW, MAIN 
SPROUT, BIRTH 


GGNORD 
GROW 


GROW 


Table 5—Parameters 


Source of values 


Name Definition Units or relationships Subroutines 
AAA(i,3) Array of species names KILL, GROW, INPUT 
A\(i) Intercept for stemwood kg Gholz and others (1979) KILL, INPUT, GROW, OUTPUT 
A2(7) Intercept for the foliage kg Gholz and others (1979) KILL, INPUT, GROW, 
AL Available light (See GROW) GROW 
AREAL Total leaf area m2 /m2 (See OUTPUT) OUTPUT 
AToT Total number of trees in plot (See OUTPUT) OUTPUT 
BAR(i)1 Total aboveground wood t/ha Gholz and others (1979) OUTPUT 
and foliage biomass 
BB(i) Slope for stemwood biomass Gholz and others (1979) INPUT, KILL, GROW, OUTPUT 
BC(7) Slope for total foliage biomass Gholz and others (1979) OUTPUT, KILL, GROW, INPUT 
equa tion 
By (i) Coefficient relating tree Reed and Clark (1979) INPUT 
volume growth to leaf biomass 
Bo(i)} Height parameter Botkin and others (1972) BIRTH, INPUT, KILL, GROW 
B3(7)1 Height parameter (calculated from 
DBHMX and HMAX) 
(see GROW) 
Ba(i) Exponential relating tree Reed and Clark (1979) INPUT 
yolume growth to leaf biomass 
Cli) Conversion from leaf mass 22/kg Gholz and others (1979) KILL, GROW, INPUT 
to area 
DBHMX Maximum diameter at breast cm Franklin and Dyrness (1973) KILL, GROW, INPUT 
height 
DIMAX (i) Maximum number of degree days Emmingham and Waring KILL, SELECT, GROW, 
growing days (Personal communication) INPUT 
pc! Diameter increment shading cm Modified from Botkin GROW 
and competition and others (1972) 
(see GROW) 
FBIol Total foliage biomass kg Gholz and others (1979) GROW, OUTPUT, BIRTH 
FI Fire intensity (O<FI<1) DISTRB 
HMAX (i) Maximum height of species i dm See Table 7 INPUT 
HT1 Tree height dn (See GROW) GROW 
IGEOG Geographic region INPUT, GROW, SELECT, KILL 
HT Decimeter height class for trees (See GROW) GROW 
IREG(j ,IGEOG) Indicates species availability KILL, SELECT, BIRTH, INPUT 
in geographic region IGEOG 
ISEL(i) Summarizes whether species KILL, SELECT, BIRTH, GROW 
a. can occur in a region 
ITOL(i) Indicates shade tolerance (See GROW) KILL, GROW, INPUT 
MPLANT. Randomly chosen number of BIRTH 
species to plant (0 to 8) 
NCODE( i) Species code number BIRTH, GROW, SPROUT, OUT- 
a; PUT, KILL, PLTREE 
NTOT Total number of trees in plot KILL, SPROUT, GROW, OUT- 
PUT, BIRTH 
P Probability a tree survives (See DISTRB) DISTRB 
a disturbance (wind or fire) 
PR Probability a tree species (See KILL) KILL 
survives 1 year of slow 
growth 
WSPRT Randomly selected number of BIRTH, SPROUT 
trees to sprout 
ssiol Total foliage biomass and stem t/ha Gholz and others (1979) GROW 
wood biomass (without bark) (See GROW) 
except for alder 
SLAR Shading leaf area, leaf above w2/m2 (See GROK) GROW 
that of size (dm) of tree 
SPNAM( 1) Alphanumeric species names 
= being considered 
SPRTMN (i) Minimum d.b.h. for sprouting to cm KILL, GROW, INPUT 
occur 
SPRTMX(i) Maximum d.b.h. for sprouting cn DBHMX except for alder KILL, GROW, INPUT 
SPRTND (7) Average number of sprouts on SPROUT, KILL, GROW, INPUT 
death of main stem 
SUMLA(IHT) Sum of leaf area of trees w2/m2 GROW 


in dm size class IHT 
(latter of all dm classes) 


WAIN (1) Minimum value of water stress bars SELECT, GROW, KILL, INPUT 
WATDEX Moisture index {See GROW) GROW 
Z(K) Normal distribution (See GGNORD) GROW, GGNORD 


Icalculated parameters. 


Table 6—Internal parameters 


1 


Name Definition Subroutine 
KSPRT(i) = -1 if sprouting occurs KILL, SPROUT, INIT 
= 1 if no sprouts 
KLAST =KTIMES, current iteration MAIN 
KYR Current year (0 KYR NYEAR) MAIN 
NCT Keeps track of number of years MAIN 
in the main program 
KYR1 KYR = 1 MAIN, OUTPUT, 
BIRTH 
NG Number of trees not growing KILL, GROW 
Nl NSPEC + 1 OUTPUT 
N2 NSPEC + 2 OUTPUT, NSPC 
NSPC Species to sprout iteration SPROUT 
NUM Species number (1 = NUM = NSPEC) INPUT 
NSP Randomly selected species to plant BIRTH 
NEW Possible species to sprout KILL, SPROUT 
NEWTR Species of possible new trees KILL, GROW 
NOGRO( i ) Equals 0 if tree is not growing INIT, KILL, GROW 
RANDOM Random number 
RAT Random number which is fixed for BIRTH 
year 
ZNYR KLAST MAIN 


Ithose that keep track of processes in the FORTRAN program itself. 


Table 7—Parameter values for CLIMACS 


2 
Growth Parameters~ 


Species SPRTMN SPRTMX! —— KTIME Al A2 BB BC c MINS 
Abies amabilis 0 0 9999 -3.5057 -4.5490 2.5744 2.1926 15.2000 18 
Abies grandis 0 0 9999 -3.7389 -3.4660 2.6825 1.9278 13.1000 18 
Abies lasiocarpa 0 0 9999 -3.7389 -3.4660 2.6825 1.9278 13.2000 10 
Abies procera 0 0 9999 3.7158 -4.8730 2.7592 2.1683 13.2000 18 
Acer macrophy]1um 10 50 9999 -3.4930 -3.7650 2.7230 1.6170 26.1000 18 
Alnus rupra 10 75 9999 0 0 0 0 28.2000 8 
Arbutus menziesii 10 120 9999 -3.7080 -3.1230 2.6580 1.6930 12.6000 30 
Castanopsis chrysophylla 10 90 9999 -3.7080 -3.1230 2.6580 1.6930 12.6000 20 
Chamaecyparis TAR eGTElC 0 0 9999 -2.0927 -2.6170 2.1863 1.7824 17.7000 10 
Libocedrus decurrens 0 0 9999 -2.0927 -2.6170 2.1863 1.7824 17.7000 27 
Picea engelmannii 0 0 9999 -3.7389 -3.4660 2.6825 1.9278 17.4000 8 
Picea sitchensis 0 0 9999 -3.7389 -3.4660 2.6825 1.9278 17.4000 5 
Pinus contorta 0 0 9999 -2.9849 -3.6190 2.4287 1.8362 17.4000 15 
Pinus Tambertiana 0 0 9999 -3.9840 -4.0230 2.6667 2.0327 17.4000 25 
Pinus monticola 0 0 9999 -4.2847 -3.9740 2.7180 2.0039 17.4000 10 
Pinus ponderosa 0 0 9999 -4.4907 -4.2610 2.7587 2.0967 17.4000 20 
Pseudotsuga menziesii 0 0 9999 -3.0396 -2.8460 2.5951 1.7009 17.4000 25 
Quercus garryana 10 90 9999 -3.4930 -3.1230 2.7230 1.6930 12.6000 25 
Thuja plicata 0 0 9999 -2.0927 -2.6170 2.1863 1.7824 17.7000 10 

Suga heterophylla 0 0 9999 -2.1720 -4.1300 2.2570 2.1280 20.4000 18 
Tsuya mertensiana 0 0 9999 -4.8164 -3.8170 2.9308 1.9756 20.4000 10 


ISPRTMX = DBHMX except for Alnus rubra which does not sprout from large trees. 


2Gholz and others (1979). 


3Zobel and others (1976), Waring (1969), and Waring, personal communication. 


Table 8—Parameter values for moisture index 


Degree Days Height Parameters 

SPECIES Max. Min. B3 B2 ITOL AGEMX SPRT ND 
Abies amabilis 80 20 0.22 85.60 5 450 0 
Abies grandis 100 50 -21 76.05 3 350 0 
Abies lasiocarpa 60 20 -78 95.18 3 350 0 
Abies procera 70 30 +19 68.29 2 450 0 
Acer macrophy11lum 115 85 ~55 55.48 3 250 3 
Alnus rubra 90 45 54 81.68 2 150 3 
Arbutus menziesii 130 75 -16 38.35 2 200 3 
Castanopsis chrysophylla 130 60 14 74.73 2 150 3 
Chamaecyparis nootkatensis 60 20 628} 69.29 4 1500 0 
Libocedrus decurrens 130 60 -30 72.72 4 750 0 
Picea engelmannii 60 20 -46 92.94 2 650 0 
Picea sitchensis 90 60 -14 64.03 3 450 0 
Pinus contorta 70 20 dey 77.63 1 200 0 
Pinus lambertiana 130 70 34 85.81 3 450 0 
Pinus monticola 70 40 -16 63.79 3 450 0 
Pinus ponderosa 130 100 -38 93.81 1 650 0 
Pseudotsuga menziesii 130 30 -23 80.44 2 1000 0 
Quercus garryana 130 100 -29 52.51 1 350 3 
Thuja plicata 90 45 -23 64.71 4 750 0 
Tsuga heterophylla 95 35 ot 81.21 5 600 0) 
Tsuga mertensiana 55 20 -07 75.97 4 750 0 
Table 9—Parameter values for shade tolerance equation 

‘ 1 2 3 2 3 
Species IREG DBHMX Bl HMAX B4 DRANGE 
1 2 3 4 

Abies amabilis 1 1 1 1 183 15612 6096 0.9 G45) 
Abies grandis 1 1 0 0 183 12493 7096 1.0 1.31 
Abies Tasiocarpa 1 a at 1 61 7404 3040 1.0 1.25 
Abies procera 1 1 1 0 220 54435 7920 of) 0.82 
Acer macrophy1]lum 1 1 1 1 50 19044 1524 1.2 -37 
Alnus rubra 1 1 1 1 75 130 1320 2.0 .24 
Arbutus menziesii 1 1 0 1 120 5047 2438 ial 1.38 
Castanopsis chrysophylla 1 0 0 0 90 15192 3500 eal 1.38 
Chamaecyparis nootkatensis 1 1 1 1 150 3135 5334 ial 1.73 
Libocedrus decurrens 1 0 0 0 120 2243 4500 1.0 1.73 
Picea enyelimannii 1 1 1 10) 102 3313 4877 9 1.31 
Picea sitchensis 0 0 1 1 230 12391 7500 1.0 1.31 
Pinus contorta 1 1 1 al 75 6650 3048 1.0 -78 
Pinus Tambertiana 1 0 0 0 125 12233 5500 9 LS 
Pinus monticola 1 1 1 1 205 15645 6675 9 1.04 
Pinus ponderosa 1 0 0 0 125 14377 6000 9 1.14 
Pseudotsuga menziesii 1 1 1 it 220 13193 9000 ibGal 1.00 
Quercus garryana 1 1 0 0 90 4956 2500 taal -74 
Thuja plicata 1 1 1 1 300 3186 6000 1.0 1.73 
Tsuga heterophylla 1 1 1 1 122 17645 6858 “9 1.47 
Tsuga mertenstana 1 1 1 1 102 9098 4572 1.0 1.11 


lyalues from Fowles (1965). 


2values from Franklin and Dyrness (1973), Fowles (1965), and Pardo (1973). 


3Adjusted forin values yiven in Reed and Clark (1971) so that a given 03H increment is reasonable. 


Subroutine BIRTH The subroutine BIRTH is based on the FORET model of Shugart and West (1977) 
(figure 3 presents a flow chart for the routine). Species available for regeneration 
are selected depending on the existing water stress, temperature range, and geo- 


graphic region of the plot. 


Select species 
(randomly) 
Establish 0 to 8 
trees, set ages to 0 
Calculate d.b.h. 
of new trees 


Biomass 
below minimum 
stocking? 


Set minimum 
d.b.h. 


Set NPLANT to 
determine number 
of cycles through 
BIRTH 


| Calculate total 
foliage biomass 


Check for species 
that can germinate 


Water 
stress, 
temperature 
range and geographic 
region appropriate 
for species? 


Yes 


No 


Last 
loop through 
NPLANT?. 


Set new 
tree ages 
to 1 


Can 
sprouting 
occur? 


Call 
SPROUT 
No 


Determine number of 
individuals to plant 


Figure 3.—Flow chart for subroutine BIRTH. 


10 


Regeneration depends on the projected leaf area of any existing trees. During 

the first year any available species can regenerate. Thereafter if the projected leaf 
area is greater than 3 m?/m? (square meter of leaf area per square meter of 
ground), Pinus contorta Dougl. ex Loud., P. ponderosa Dougl. ex Laws., and 
Quercus garryana Dougl. ex Hook. cannot regenerate. If the projected leaf area is 
greater than 6 m?/m?, then Abies procera Rehd., Alnus rubra Bong., Arbutus men- 
Ziesii Pursh, Castanopsis chrysophylla (Dougl.) A. DC., Picea engelmannii Parry 
ex Engelm., and Pseudotsuga menziesii cannot appear in the plot. For a projected 
leaf area exceeding 8 m?/m? Abies grandis (Dougl. ex D. Don) Lindl.), A. /asio- 
carpa (Hook.) Nutt., Acer macrophylilum Pursh, Picea sitchensis (Bong.) Carr., 
Pinus lambertiana Dougl., and P. monticola Dougl. ex D. Don cannot appear. If 
the projected leaf area exceeds 10 m?/m?, Libocedrus decurrens Torr., Chamae- 
cyparis nootkatensis (D. Don) Spach, Thuja plicata Donn ex D. Don, and Tsuga 
mertensiana (Bong.) Carr. cannot regenerate. At leaf areas exceeding 10 m?/m? 
only the most shade tolerant species, Abies amabilis Doug!. ex Forbes and Tsuga 
heterophylla (Raf.) Sarg., can regenerate. The trees appear in the plot when they 
are between 10 and 15 cm d.b.h. This limit on size is necessary because num- 
bers of smaller seedlings become so large that computer time is prohibitive, and 
the abundance of very small seedlings may not reflect the composition of an 
understory that could replace large-canopy trees. 


The projected leaf area (XLAI) (in m?/m?) is based on foliage biomass (FB10) of all 
trees on the plot. The equations for the foliage biomass of trees less than 50 cm 
(19.7 in) in diameter are from Binkely (1983) for A/nus rubra and from Gholz and 
others (1979) for all other species. Because the data range of Gholz’ equations is 
exceeded for large trees (DBH > 50) (fig. 4), the equations for the foliage bio- 


320 
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Figure 4.—Foliage biomass as related to diameter using the 

equations of Gholz and others (1979). A warning is given by 

Gholz and others (1979) that “maximum and minimum values 

from the data set...do not necessarily mean that the equation 

provides reasonable numbers at these extremes.” Because of 

problems in exceeding the data ranges, the foliage biomass- 

d.b.h. relationship was modified using figure 5, as described 

in the text. 11 


mass of large trees are calculated indirectly. The sapwood-diameter equation for 
Pseudotsuga menziesii has been obtained by regression analysis.~ Foliage bio- 
mass can then be calculated from the foliage biomass-sapwood area ratio (Grier 
and Waring 1974) (fig. 5). To estimate the foliage biomass for large trees of other 
species, a multiplier is calculated from the ratio of the foliage biomass at 50 cm 
diameter of the species being considered to that of P. menziesii at 50 cm diam- 
eter. The multiplier is used to adjust the foliage biomass equations from Gholz 
and others (1979) for trees larger than 50 cm. 


‘Personal communication, R. H. Waring, School of Forestry, 
Oregon State University, Corvallis, Oregon. 
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Figure 5.—Foliage biomass as related to diameter for Pseudo- 
tsuga menziesii. 


The equations are as follows: 


exp[Ao(/) + BC(/) In DBH(k) ] DBH(k) <50, i; # 6 
FBIO = { 3.2 + 1.89 In (DBH(k)), DBH(k) =50, /= 6 
[256 * In (DBH(k)) — 639] *2.5*DRANGE(/)/C(16), DBH(k) >50 ; 


where: 


i = 1,....NSPEC (NSPEC = total number of species) 
(i = 6 = Alnus rubra), 


ij = 1,...NTOT (NTOT = total number of trees), 


A, (i 
| = parameters from Gholz and others (1979) (see table 7), 
i 
DBH(k) = diameter at breast height of tree k (cm), 
DRANGE(/) = the ratio of the foliage biomass of species / to that of 
P. menziesii at 50 cm, and 
C(16) = conversion from kg to m? for P. menziesii. 


To calculate the leaf area for each tree the foliage biomass is multiplied by the 
factor C(/) to convert from kg to m? and is divided by 2.0 for deciduous trees and 
2.5 for conifers to convert total leaf area to projected leaf area. The projected 
leaf areas are then summed for all trees and divided by the plot size (2000 m?) to 
obtain a value for square meter of leaf area per square meter of ground surface: 


aay = eGo eBlomy 
k=1 2000 ; 
2.0 if / is a deciduous tree and 
where: P(/) = Aas : 
2.5 if / is a coniferous tree . 


The projected leaf area of the plot is tested against the maximum leaf area at 
which available species can regenerate. If the projected leaf area is less than the 
species maximum, a random number of saplings (from 1 to 8) enters the plot. 
Each introduced sapling has a randomly selected diameter of 10 to 15 cm (3.9 to 
5.9 in). 


Sprouting may occur depending on tree size (subroutine SPROUT). Species that 
reproduce vegetatively in this manner are Acer macrophyllum, Alnus rubra, 
Castanopsis chrysophylla, Quercus garryana, and Arbutus menziesii. 


BIRTH continues until the projected leaf area exceeds 1.0 m?/m? resulting in full 
stocking the first year of the model run. 


2/Leaf area conversions from personal communication, 
R. H. Waring, School of Forestry, Oregon State University, 
Corvallis, Oregon, based on Grier and Waring (1974). 


Subroutine GROW 


Optimum Growth Per Year 


Decreases in Growth 
Rate Related to Size 


Subroutine GROW calculates diameter increment for each tree (fig. 6). The annual 
diameter increment function is similar to that used by Botkin and others (1972), 
Shugart and West (1977), and Mielke and others (1978). CLIMACS is unique in 
that it incorporates a specific function for the effect of soil moisture upon growth. 
Diameter growth takes into account the following conditions: 


DaARONMH 


. Growth per year is dependent upon existing foliage biomass 
. Growth decreases as tree size increases 
. Temperature affects photosynthetic rate 


Available moisture affects growth 
Competition influences growth 
Shading and shade tolerance affect growth. 


The growth subroutine calculates the maximum growth of a species and modifies 
that optimal growth by size and site considerations. There is strong evidence that 
diameter growth is related to existing foliage biomass (Grier and Logan 1977). 
CLIMACS uses a nonlinear relationship between leaf biomass and growth: to cal- 
culate optimal diameter increment B: (i) FBIO®4". The parameters for this 
relationship are adjusted from those of Reed and Clark (1979) to give reason- 
able diameter increments for trees 50 cm (19.7 in) in diameter. 


As a tree increases in size, an increasing portion of its photosynthetic output is 
used in respiration and not in producing additional tissue. The phenomenon of 
decreasing growth in large trees is represented by the following calculation 
from Botkin and others (1972):~ 


_ DBH() HT (/) 


’ 


DBHMX(/) HMAX(/) 


for species / and tree /. 


As a tree becomes very large, growth approaches zero. Maximum diameter and 
height were chosen for each species to represent large but not the largest trees 
of a species to avoid record-sized stands. 


2’Symbols are defined in table 5 and parameter values given 
in table 7. Because growth is conditioned by the maximum 
diameter and height for the species, the maximum values 
determine the upper limit for tree growth. When using the 
data for record trees (Fowles 1965, Pardo 1973) the model 
always grows Site 1 trees. The current values for DBHMX(/) 
and HMAX(/) are for large but not the largest trees of a spe- 
cies (from Franklin and Dyrness 1973). 


Calculate foliage biomass 
of this tree 


Calculate maximum, open- 
growth diameter increment 
F (DBHMX,FBIO,DBH) 


Set minimum growth 
limit, set available 
light constant 


Calculate reduced diameter 
increment for tolerant, inter- 
mediate, and intolerant species 
DINC F (DNC,AL) 


Reset supression flags 
to O, INITIALIZE leaf 
area sum and biomass sum 


Calculate total stand, 
stem and foliage biomass 
Calculate Height leaf area Set suppression flag 
profile by decimeter intervals 


Next tree 


Increment 
less than minimum Increment 
set? DBH 


Calculate growth moisture 
index (dry, intermediate, 
or wet) 


Calculate tree height and 
shading of trees taller 


Figure 6.—Flow chart for subroutine GROW. 


Table 10—Regression values for height as compared to diameter where 
HT = 137 + B2 DBH + Bs DBH* * 


Stressed or Bo B3 Re N 
Species nonstressed!/ D.b.h. range 
Centimeters 
Pseudotsuga menziesii Stressed 72.34 -0.22 0.98 55 9.8 - 196.4 
Nonstressed 80.44 -.23 -96 65 22.3 - 200.4 
Tsuga heterophylla Stressed Uo sar2// -98 38 4.9 - 143.4 
Nonstressed S121 = 5277, -97 39 6.8 - 130.0 
Abies procera Stressed 68.29 -.19 -98 24 48.1 - 127.3 
Abies amabilis Nonstressed 85.60 -.55 -98 19 16.7 - 62.0 
Tsuga mertensiana Nonstressed 55.29 -.07 -98 23 32.2 - 84.2 
Thuja plicata Nonstressed 64.71 -.23 -98 38 6.2 - 157.9 


*Data is from research stands in the H. J. Andrews Experimental Forest. 


1/stressed trees are growing in a habitat where PMS < 15 or TGI > 60. 


The height-diameter relationship used for some species is based on that used 
by Botkin and others (1972) and Mielke and others (1978). The parabolic equa- 
tion follows from a discussion in Ker and Smith (1955). Height (HT) is a func- 
tion of the diameter at breast height (DBH): 


HT(j) = 137 + Bo(i) DBH(/) — Ba(’) DBH(/)°. 


The constants Ba(/) and Bsa(/) are chosen so that height is at a maximum and the 
rate of change of height with respect to diameter is zero (dHT(/)/dDBH(/) — 0) 
when DBH(/) is at a maximum: 


Bo(i) = AIMAX() 187) and 
DBHMX(/) 


we a IMAX(i) aa 13% 
BSW) Steen grey era ee 
DBHMX(/) 
The rate of change in height decreases as the diameter increases. For Tsuga 
heterophylla, T. mertensiana, Thuja plicata, Pseudotsuga menziesii, Abies ama- 
bilis, and A. procera height was estimated directly from diameter using equa- 
tions developed by regression analyses of diameter-height data from the 
H.J. Andrews Experimental Forest (Willamette National Forest, Oregon) (table 10). 
Sufficient data were not available to treat all species in this manner. Because 
height growth is related to habitat and Tsuga heterophylla and P. menziesii were 
sampled in a variety of habitats, regressions for these two species were also per- 
formed for temperature- or moisture-stressed stands (TGI > 60 or PMS < 15) 
and nonstressed stands. In stressed stands the height growth is less—particu- 
larly for large trees (fig. 7). 
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Figure 7.—Tree height as related to diameter for Tsuga hetero- 
phylla and Pseudotsuga menziesii for stressed and nonstressed 
stands from regression analysis. 


Temperature 


Effects of Moisture 


The effect of temperature on the photosynthetic rate is modeled with a parabolic 
equation that ranges from 0 to 1 (following Botkin and others 1972): 


4(TGI — DMIN(/)) (DMAX(i) — TGI) 


TEGD = ; 
(DMAX(i) — DMIN(i) ) 


A value of 0 occurs for TGI less than the minimum degree days (DMIN(/)) or greater 
than the maximum degree days (DMAX(/)). Thus, this function reflects the tem- 
perature range in which a species can grow. The values for DMIN(i/) and DMAX(/) 
are from Emmingham and Waring.” 


The inclusion of the moisture index (WATDEX) in the equation for diameter incre- 
ment is necessary due to the major influence of available moisture on the growth 
of tree species in the Pacific Northwest (Franklin and Waring 1980). In this model 
the effect of available moisture on growth is represented by a beta function (fol- 
lowing Reed and Clark 1979): 


Ote)), (CS0Y 


, for x€(a,c) 
WATDEX = B(a,b,c,x) =$ (8-4) (c—b) 
0, otherwise ; 
where: x = PMS, 
c = WMIN(/) = the minimum value of water stress for a 


species (—bars), 
a = the negative x-intercept, 
b = the optimum PMS for growth (—bars), and 
v = (c—b)/(b—a). 


See table 11 for a list of the parameter values. 
“Personal communication, W. Emmingham and R. H. Waring, 


School of Forestry, Oregon State University, Corvallis, 
Oregon. 


Table 11—Parameter values for moisture index' 
Moisture a b Species 


Dry -10 5 Abies grandis, Arbutus menziesii, 
Acer macrophyllum, Castanopsis 
chrysophylla, Libocedrus decurrens, 
Pinus contorta, P. Tambertiana, P. 
ponderosa, Pseudotsuga menziesii, 
Quercus garryana 


Intermediate -8 1 Abies amabilis, A. lasiocarpa, A. 


procera, Picea engelmannii, Pinus 
monticola, Thuja plicata, Tsuga 
heterophylla, T. mertensiana 


Moist -2 0 Alnus rubra, Chamaecyparis 
nootkatensis, Picea sitchensis 


1The parameters, a and b, are explained in the text. 


According to the model, dry site species have a high value of WATDEX for most 
values of PMS (see fig. 8). Intermediate site species, such as Tsuga heterophylla, 
grow only at lower PMS values, and the wet site species, such as Picea sitchensis, 
have a positive value of WATDEX only at very low values of PMS. 


‘A 
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® Pseudotsuga menziesii 
4 Abies amabilis 
@ Alnus rubra 


Proportion of optimum growth 


0 30.00 
Plant moisture stress (-bars) 


Figure 8.—The effect of plant moisture stress upon diameter 
increment. ; Fimant Gina, ; 
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Competition for Nutrients 


Shade Tolerance 
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The crude expression for rooting space competition, which reflects competition 
for nutrients, is from JABOWA (Botkin and others 1972) and is dependent upon 
the total aboveground biomass (TBIO) and the maximum biomass a stand can 
support (SOILQ): 1 — TBIO/SOILQ. The total aboveground biomass (TBIO) is the 
sum of the foliage biomass and stemwood biomass. The equations for stemwood 
biomass are from Gholz and others (1979): 


0.02 + 2.09 DBH(j) — .0015 DBH(j)*, i = 6 (Alnus rubra) 


SBIO -| f 
exp(Ai(/) + BB(/) In DBH(/)) , i~6; 
where A;(/) and BB(/) are parameter values from Gholz and others (1972) (table 7). 
The data range for the Gholz equations for stemwood biomass does not include 
large trees. To assess the impact of this problem stemwood biomass was plotted 
against values from yield tables (MacLean and Berger 1976) (fig. 9). Modeled and 
measured biomass values for Calocedrus decurrens and Pinus ponderosa com- 
pare well, but the predicted curve for Pseudotsuga menziesii stemwood biomass Is 
greater than measured biomass at large diameters. As there are no data available, 
evaluation of this discrepancy is not possible. 


The maximum biomass for a site (GOILQ) is obtained from recorded values of 
naivta! stands similar to the habitat for which the model is being run. As TBIO 
approaches SOILQ from below, (1-TBIO/SOILQ) approaches 0, and diameter 
increment decreases. When TBIO exceeds SOILQ, (1-TBIO/SOILQ) becomes 
negative; no diameter increase occurs. The value for SOILQ should be chosen 
carefully because evaluation of the model output suggests that the stand biomass 
is strongly influenced by this value. 


Diameter growth is also related to the light a tree receives and its shade tolerance. 
The regulation of photosynthetic rate by available light is expressed as: 


1 — exp(—t(AL — f2)) ; 


where ti and fe are species specific parameters (see table 12). The available light 

a tree receives (AL) is based on the total foliage area of all trees in the plot that 
are taller than the tree being considered (fig. 10 illustrates how available light 
changes over time). The canopy is broken into tree height intervals by decimeters. 
The sum of the leaf areas of all trees in height classes greater than the one being 
considered is calculated and stored as the shading leaf area (SLAR). Available 
light is then found by using an equation modified from Botkin and others (1972) and 
Mielke and others (1978): 


AL = PHI * exp(—0.26 SLAR) ; 


where PHI = 1.0 is the annual insolation (in appropriate units). 
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Figure 9.—The relationship of stemwood biomass to diameter 
for (A) Calocedrus decurrens, (B) Pseudotsuga menziesii, and 
(C) Pinus ponderosa. The curves are from Gholz and others 
(1979), and the line crossing the curve indicates the end of 
the data range. Circles are from volume tables in MacLean 
and Berger (1976) with solid circles being inside and open 
circles outside the range of sampled trees. 
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Table 12—Parameter values for shade tolerance equation 


20.467 


32.846 


40.0 


te 


-025 


-025 


-025 


-025 


Shade tolerance class based on 
Minore (1979) 


Quercus garryana, Pinus ponderosa, P. 
contorta 


Acer macrophyllum, Alnus rubra, 
Castanopsis chrysophylla, Arbutus 
menziesii, Picea engelmanni, Abies 
procera 


Abies grandis, A. lasiocarpa, Pinus 


Tambertiana, P. monticola 

P. sitchensis, Pseudotsuga menziesii, 
Libocedrus decurrens, Chamaecyparis 
nootkatensis, Thuja plicata, T. 
mertensiana 


Abies amabilis, Tsuga heterophylla 


ITthe relative diameter growth is equal to (1-exp (-t, (AL-to))) for AL 
being the available light derived from the projected leaf area (m2/m2). 


50 Year 1 Light 


Height (m) 
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y Pseudotsuga Menziesii 

x x Castanopsis chrysophylla 
AW a/ Pinus ponderosa 
----Calocedrus decurrens 


Foliage biomass (kg) 


Figure 10.—The relationship of height to foliage biomass from 
a simulation model. The size of the arrow representing light is 
meant to depict the amount of light reaching the canopy at 
each level. 
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Diameter Increment 
Equation 


24 


CLIMACS includes five curves that translate available light into photosynthetic 
capacity based on shade tolerance (fig. 11). Pinus ponderosa is the least tolerant 
of the Pacific Northwest species and cannot grow under a leaf area greater than 
5 m2/m?2 (Minore 1979; Waring). Hardwoods, intolerant conifers, and tolerant 
conifers, respectively, represent increasing levels of shade tolerance (Minore 1979). 
In the model the greatest difference in these categories occurs at a projected leaf 
area of 8 m?/m?. 


The diameter increment equation is a multiplicative function of six factors: 


DINC -( a | DBH(j) B:(i) FBIO®4"” (1-TBIO/SOILQ)TEGD 
DBHMX(i) HMAX(i) 274 + 3B2(i) DBH(/)* —4B3(/) DBH(/)* 


x (WATDEX) X [1 — exp(-a1(AL-az2) )] ; 


where: / = 1,...,.NTOT = the tree number and 
i = 1,....NSPEC = the tree species code. 


This assumes that optimum growth is modified by the effects of temperature, 
moisture availability, competition, and shade tolerance, all acting independently. 
A multiplicative function is the most common relationship used (Botkin and 
others 1972, Reed and Clark 1979, Shugart and West 1977) but can tend to under- 
estimate growth if many factors are incorporated (Swartzman and Bentley 1979). 


2/Personal communication, R. H. Waring, School of Forestry, 
Oregon State University, Corvallis, Oregon. 
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Figure 11.—The effect of shade tolerance on diameter growth. 


Subroutine KILL Much of the existing literature on mortality predicts death of a tree based on its 
crown quality or vigor class (Graham 1980, Hamilton and Edwards 1976, Staebler 
1953). Because we cannot identify crown features in the current model structure, 
individual tree mortality is based on diameter, the maximum diameter for the spe- 
cies, and the life history of the species (see fig. 12). A tree is subjected to KILL if 


Short-lived Long-lived Long-tived 
early seral early seral mid-searal mid-seral 
witha 
increasing 
mortality 


Yes 
Check for 
spouting 


Update diameters, 
ages and species 
codes 


ee 


Yes 


Figure 12.—Flow chart for subroutine KILL. 
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Table 13—Equations for slow growth-related mortality 


Successional status Conditions Probability of surviving 1 year Species 

Long-lived DBH < 0.2DBHMX 0.8866 + 0.1029 DBH-1 Pinus ponderosa 

early seral 0.2-1/DBHMX DBHMX Pseudotsuga menziesii 
DBH > 0.2DBHMX 0.9895 

Short-lived 0.628 Acer macrophy1lum 

early seral Alnus rubra 


Arbutus menziesii 
Quercus garryana 


Long-lived DBH < 0.1 DBHMX 0.96496 + 0.01294 DBH-1 Chamaecyparis nootkatensis 
mid-seral -1-1/DBHMX  DBHMX Libocedrus decurrens 
huja plicata 
DBH > 0.1 DBHMX 0.97556 + 0.02339 DBH 
= DBHMX 
Long-lived DBH < 0.5DBHMX 0.8866 + 0.1112 DBH-1 Abies grandis 
mid-seral ne ie 0.5-1/DBHMX  DBHMX he lasiocarpa 
increasing mortality + procera 
Picea engelmannii 
DBH > 0.5DBHMX 1.0538 - 0.1112  DBH P. sitchensis 
= “DBHMX —t*” Pinus contorta 
P. Tambertiana 
P. monticola 
Late seral DBH < 0.1 DBHMX 0.9956 + 0.99335 DBH-1 Abies amabilis 
0.1-1/DBHMX DBHMX suga heterophylla 
J. mertensiana 


DBH > 0.1 DBHMX 1.0053 - 0.0635 DBH 
isi ~DBHMX 


Compiled from McArdle and others (1949) and from unpublished data from a chronosequence of stands in the 
Cascade Range in Oregon and Washington (Ted Thomas, USFS, Forestry Sciences Laboratory, Corvallis). 


diameter growth is less than 1 mm (0.04 in) in any year. We classified tree species 
into five groups based on mortality pattern using data from a chronosequence of 
stands in the Cascade Range (table 13). Short-lived, early seral species have a 
constant survival probability of 0.628 for 1 year. For all other successional classes, 
mortality changes with the tree size (fig. 13). Long-lived, early seral species have 
increasing survival until 20 percent of the maximum diameter is reached, at which 
point the probability of dying becomes constant. Late seral species have an in- 
creasing probability of survival until 10 percent of the maximum diameter is at- 
tained; then the probability of surviving tapers off slowly. For long-lived, mid-seral 
species the probability of survival slowly increases as size increases. Long-lived, 
mid-seral species with increasing mortality have a peak of surviving at 50 percent 
of the maximum diameter; the probability of death increases on either side of this 
asymptote. 
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Figure 13.—Probability of tree survival 1 year after slow 
growth (less than 1 mm [0.04 in] of diameter increment) for four 
successional classes. A fifth class, short-lived early-seral 
species, has a constant survival probability of 0.628. 


The subroutine DISTRB models natural and human disturbances to the eco- 
system. Fires are the most common natural disturbances to forests of the Pacific 
Northwest (Hemstrom 1979), and windthrows can cause localized disturbances. 
Clearcuttings are simulated on user-specified rotation bases with Pseudotsuga 
menziesii planted after the cut (other species could be planted). Each of these 
disturbances can be introduced at user-specified intervals and intensities (see 
subroutine INPUT, Appendix 1). 


Because the BIRTH subroutine results in small trees (that is, the seedling stage is 
skipped), model year 1 is approximately stand year 30 (assuming a 10- to 15-cm 
(3.9- to 5.9-in) tree is found on a 30-year-old-stand). A disturbance should be 
thought of as occurring approximately 30 years in the past. Thus a disturbance 
interval of 100 years results in the disturbance occurring in model years 100, 200, 
300, etc., or in stand years 130, 260, 390, etc. 


Following an episodic disturbance, the species, the diameter, and the age of all 


trees killed can be printed out in tabular form and can be saved on a tape for 
further use (refer to Appendix 1, subroutine INPUT). 
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Fire is modeled as an event that causes a high proportion of death for young trees 
dependent upon the species’ tolerance to fire (fig. 14). Pseudotsuga menziesii, 
Pinus ponderosa, P. lambertiana, Picea sitchensis, Abies grandis, A. lasiocarpa, 

A. procera, and Quercus garryana are the species most resistant to fires (Minore 
1979; Franklin®); and Abies amabilis, Tsuga heterophylla, and T. mertensiana 
usually don’t survive fire The probability of a tree surviving a fire (Pr) is related 
to the intensity of the fire Fl (0 = FI < 1): 


exp(—0.00255*DBH(/) )*Fl i = 2, 3, 4, 12,14, 15, 16, 17 
Pr= {4 exp(-0.00053*DBH(/))*Fl, / = 5,6, 7, 8, 9, 10, 11, 13, 18, 21 
0 AO 2ZOn 


Fire occurs at user-specified intervals. The present model probably works for 
severe fires that result in stand replacement but is inappropriate for the more fre- 
quent and less intense fires in forests dominated by Quercus garryana or Pinus 
ponderosa. 


The equations for the probability of death of a tree by a catastrophic wind are 
based on data from Tar Creek in the Clackamas River drainage (see footnote 6). 
Individual Pseudotsuga menziesii, Tsuga heterophylla, and Thuja plicata that died 
as a result of wind between 1948 and 1958 were recorded. All species have an 
increasing probability of death from wind as size increases (fig. 15). Although 
large trees have a large surface area upon which the wind can act, smaller trees 
are frequently in the path of falling large trees. The equations in the model are 
dependent on size and are conditioned by species and size: 


P= (1.9 -exp(—0.00053*DBH(/)))*Fl, # =1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 


(1.6 —exp(—0.00255*DBH(/)))*Fl, # =11, 19, 20 

UW, Wek UE qe, Us, az, lsh, A 
A clearcutting is modeled as an event that occurs at user-specified rotation times 
and removes all of the trees. After all trees are removed from a plot, 198 trees 
per 0.20-ha (0.08-acre) plot of 10- to 15-cm (3.9- to 5.9-in) d.b.h. Pseudotsuga 
menziesii are introduced to mimic the normal planting prescription of 1,000 trees 
per hectare. 


© Personal communication, J. F. Franklin, Pacific Northwest 
Forest and Range Experiment Station, USDA Forest Service, 
Corvallis, Oregon. 


“Personal communication, J. K. Agee, National Park Service 
Cooperative Park Studies unit, College of Forest Resources, 
University of Washington, Seattle, Washington. 
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Figure 14.—Probability of death of a tree following a fire for 
fire-tolerant and fire-intolerant species. Abies amabilis, Tsuga 
heterophylla, and T. mertensiana do not survive fires in 
the model. The intensity of the fire can be adjusted by the 
parameter Fl. The most severe fire is depicted (Fl = 1). 
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1.0 Wind tolerant 


P = [1.6 — exp (—0.00255 DBH)] FI 


1.0 Wind intolerant 


Probability of dying 


P = [1.9 — exp (—0.00053 DBH)] FI 
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Diameter (cm) 
Figure 15.—Probability of death after a windstorm. Species in 
each tolerance class are listed in the text. The most severe 
storm is illustrated (Fl = 1). 


Biological Support 
of the Model 


Biological support of the model is based on three criteria: 


1. Scientific rationale of equations and parameter values, 
2. Applicability of this model paradigm to at least five forest systems, and 
3. Model output compared to data. 


CLIMACS has been verified through evaluation of the model equations and 
parameter values as explained in this text. All relationships and values are based 
on current knowledge of the tree species of the Pacific Northwest. For large trees, 
little information exists on the height, stemwood biomass, and foliage biomass 
relationships to diameter. No studies have documented the change in temperature 
growth index, plant moisture stress, or leaf area as a stand develops. 


The FORTRAN code for CLIMACS is based on models that have been applied to 
other regions. JABOWA (Botkin and others 1972) simulates forest succession in 
the northeast United States. Snugart and West (1977) used FORET, the predeces- 
sor of CLIMACS, to model the chestnut blight in the Appalachian forest of Ten- 
nessee. The FORET model has also been used to project the influence of air 
pollution on stand development (West and others 1980) and to successfully simu- 
late forest development for a mixed pine-hardwood forest in Arkansas (Mielke and 
others 1978). Other models derived from the core paradigm have been useful in 
predicting forest changes in eucalyptus forests of Australia (Shugart and 

Noble 1981) and in a tropical forest (Doyle and others 1979). Weinstein (1982) 
examined the diameter increment of a version of FORET and found that the aver- 
age increment overestimated the diameter changes of trees remeasured after a 
13-year interval. This discrepancy was corrected when moisture stress and nutri- 
ent limitation were introduced to the model. Both of these factors are in the diam- 
eter increment equation in CLIMACS. 


Simulations of forests in western Oregon compare well to existing stand structure 
(Hemstrom and Adams 1982). For xeric and mesic stands nearly 450 years old, 
CLIMACS projects a stand dominated by Pseudotsuga menziesii of a variety of 
size classes. Projected leaf area, basal area, and foliage biomass of the simulated 
stand are similar to those of the natural stand (see discussion in Hemstrom and 
Adams 1982). Differences between the natural and the simulated stand may be 
due largely to chance and the lack of disturbances in the model forest. 


Other comparisons incorporate natural disturbances in the simulation. Projections 
of 750 years of stand development on the western Olympic Peninsula of Washing- 
ton State in the aftermath of fire, windstorm, or clearcutting show the major role 
Pseudotsuga menziesii plays (Adams and others 1983). When these disturbances 
occur at various frequencies, each disturbance results in unique patterns of stand 
development (Dale and others 1984). 
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A sensitivity analysis is a process to evaluate how sensitive model output is to 
variation in parameter values. By varying parameters singly and in combination 
(for those known to be correlated), the influence of the parameter values on 
model behavior can be assessed. A sensitivity analysis has not been performed on 
CLIMACS because of a lack of time and money, and we urge that such a study 
be carried out. In the documentation process, we have noted the sensitivity of the 
model to the value for the maximum biomass for the forest (SOILQ) and to maxi- 
mum diameter for each species. The parameters affecting the diameter increment 
growth, in particular, should be examined carefully. Shugart and West (1979) 
report that FORET, a forerunner of CLIMACS, is particularly sensitive to plot size. 
Finn (1979) introduced a form factor to account for crown depth and shading of 
trees in JABOWA (Botkin and others 1972) and found little effect. Crown form 
may, nevertheless, affect the development of coniferous forests. 


As with any model, CLIMACS can be altered and improved to meet the objectives 
of the study for which it is being used. Three possible improvements to the model 
are: (1) altering the mortality curves, (2) changing the vertical structure of each 
tree, and (3) improving the competition function. 


Mortality is a largely unexplored component of successional processes. The 
model can be used to ask questions about differences in mortality among species 
by altering the mortality curves (fig. 13) and treating the survival equations as 
hypotheses to be tested by simulations. A problem with such a simulation study is 
the lack of data for comparison with model output. Data on survival rates of 
tagged trees will be available for such a study from permanent plots established 
by the USDA Forest Service. 


The vertical structure of each tree could be modified to more accurately depict 
the form of a species. This would allow analyses of the role of available habitat for 
birds and mammals as forests change over time. By correlating breeding habitat 
to vegetation structure, the probability of a breeding habitat existing under certain 
management regimes and natural successions could be assessed (for example, 
Smith and others 1979). These modifications would also allow further study of the 
impact of leaf area on growth. 


The function that models competition for nutrients in the GROW subroutine is a 
crude approximation. The model applicability would improve by deletion of the 
parameter SOILQ and incorporation of nitrogen and phosphorus cycling. The 
problem in including nutrient cycling in the model is depicting an accurate repre- 
sentation but not requiring too many driving variables. 


This research could not have been completed without the scientific advice, inter- 
est, and encouragement and unpublished data provided by Jerry Franklin. We 
appreciate technical assistance and discussions with Jim Agee, Paul Alaback, 
Bill Emmingham, Bob Harr, Mark Harmon, Tom Hinckley, Mark Klopsch, Joe 
Means, Ed Small, Gordie Swartzman, and Dick Waring. The research was sup- 
ported by contract CX-900-0-E014 from Pacific Northwest Region, the National 
Park Service, U.S. Department of the Interior. This is Publication Number 2391 of 
the Environmental Sciences Division, Oak Ridge National Laboratory, Oak Ridge, 
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Appendix 1 


Listing of 


Program CLIMACS 


(0001) 
(0002) 
(0003) 
(0004) 
(0005) 
(0006) 
(0007) 
(0008) 
(0009) 
(0010) 
(0011) 
(0012) 
(0013) 
(0014) 
(0015) 
(0016) 
(0017) 
(0018) 
(0019) 
(OG20) 
(0021) 
(0022) 
(0023) 
(C024) 
(0025) 
(0026) 
(0027) 
(0028) 
(0029) 
(0030) 
(0031) 
(0032) 
(00393) 
(0034) 
(0035) 
(0036) 
(0037) 
(0038) 
(0039) 
(0040) 
(0041) 
(0042) 
(0044) 
(0045) 
(0046) 
(0047) 
(0048) 
(0049) 
(0050) 
(0051) 
(0052) 
(0053) 
(0054) 
(0055) 
(0056) 
(0057) 
(0058) 
(0059) 
(0060) 


36 


cas 
cas 


cas 
c@s 


c@s 


FOREST 


__... FOREST SUCCESSION MODEL 


CONVERSION PROGRAMMERS: J. HOLT, M. RUST 
MILNE COMPUTER CENTER, OSU. 


All comments labeled "C@S" were done at the Center for Quantitative 
science at the University of Washington using the PRIME. To run the 
model on a CDC note these comments and make the appropriate changes. 


PROGRAM CONIFER (INPUT, OUTPUT, TAPE4, TAPES2INPUT, TAPE6=OUTPUT, TAPES 
1) 
(FILE DECLARATIONS BELOW) 


COMMON /FOREST/ NTREES(50), DBH(700), IAGE( 1000), KSPRT (50), NEWTR( 50) 
1, SUMLA( 1000), NEW(50), NCODE(700), ISEL (50) 


PRIME STORES 4 CHARACTERS PER REAL ( 1 REAL = 2 WORDS ) 
cbc STORES 10 CHARACTERS PER WORD 


COMMON /PARAM/ AAA(8, 100), DMAX(50), DMIN(50), B3(50), B2(50), ITOL(S50) 
1, AGEMX (50), 6(50), SPRTND(50), SPRTMN(50), SPRTMX(50), KTIME(50), A1 (50) 
2, A2(50), BB( 50), BC(50),C(50), WMIN( 50), IREG(S0, 4), DBHMX(50), B1(50), 
3B4(50), HMAX (50), DRANGE (50) 

COMMON /CONST/ NSPEC,. SOILQ, PMS, TGI, TGII, SD, IGEOC 

COMMON /DEAD/ NOGRO(700) 

COMMON /COUNT/ NTOT, NYEAR, IPERT. INT, IN, PMORT, FI 

COMMON /PASSX/ A(501,25),FB(50> 

COMMGN /SEED/ USEED(13) 

INTEGER USEED 


DIMENSION Z(2), SPNAM(S, 40) 
GET SOURCE DEFINITIONS OF KSREAD, K$WRIT (TO OPEN FILES) 

NOLIST /#* SYSCOM>KEYS. F 

LIST 

DATA SPNAM fp PXofoyin ys Boae) i 

DATA NCT/O/ 

INITIALIZE FILES 

CALL OPEN( 5S, ‘SUCDAT’, 6, KSREAD 

CALL OPEN( 6 ‘OUTPUT’, 6& MKSWRIT 

CALL OPEN( 7, ‘TAPE4’, 5S, KSWRIT 


CALL OPEN( & ‘TAPES’, 3S, KSWRIT 
CALL OPEN( 9, ‘TAPE9’, 5S, KSWRIT 


/* UNIT 4 RESERVED ON PRIME 


wwewewvew 


Poa READ INPUT DATA, PARAMETERS, AND VARIABLES FROM TAPES. 


(0061) 
(0062) 
(0063) 
(0064) 
(0065) 
(0066) 
(0067) 
(0068 ) 
(0069) 
(0070) 
(0071) 
(0072) 
(0073) 
(0074) 
(0075) 
(0076) 
(0077) 
(0078) 
(0079) 
(0080) 
(0081) 
(0082) 
(0083) 
(0084) 
(0085) 
(0086) 
(0087 ) 
(0088) 
(0089) 
(0090) 
(0091) 
(0092) 
(0093) 
(0094) 
(0095) 
(0096) 
(0097) 
(0098) 
(0099) 
(0100) 
(0101) 
(0102) 
(0103) 
(0104) 
(0105) 
(0106) 
(0107) 
(0108) 
(0109) 
(0110) 
(O111) 
(0112) 
(0113) 
(0114) 
(0115) 
(0116) 
(0117) 
(0118) 
(0119) 


ananana 


gaannna 


gaNnNanaANn 


10 


20 


30 


40 


42 


44 


45 


CALL INPUT 


ASSIGN FIRST TEN CHARACTERS OF SPECIES NAME TO THE ARRAY USED 
BY THE PLOTTING PROGRAM. 


DO 10 I=1, NSPEC 
CALL STRMOV( AAAC1,1), 1, 10, SPNAM(1,1), 11 ) 
CONTINUE 


IPLOT = 0 
KLAST = 1 
NYR1 = NYEAR+1 
KTIMES = 1 
N3 = NSPEC+3 
DO 20 IVA=1, NYR1 
DO 20 IVS=1,N3 
ACIVA, IVS) = O. 
CONTINUE 
CALL RANSEE 
Sse yelts PLOT SIZE IS ONE-FIFTH HECTARE. 
5 8 eG SOIL@ - THE MAXIMUM BIOMASS RECORDED FOR FORESTS IN THE AREA 
Sods THE MAXIMUM LIVE ABOVE-GROUND BIOMASS POSSIBLE IS 
eee 1000 TONS/HECTARE 
SOIL@ = 100000. 
eae TGI -— TEMPERATURE GROWTH INDEX 
Pe NSPEC — NUMBER OF SPECIES FOR SIMULATION 
CONTINUE 
CALL PLOTIN (IPLOT) 
WRITE (6,40) IPLOT, NYEAR 
FORMAT (/,15H PLOT NUMBER ,14, 24H TOTAL LENGTH OF RUN IS , 
IS, 6H YEARS////////) 
KYR = O 
CALL OUTPUT (KYR, IPLOT) 
DO 70 JJ=1,KTIMES 
CALL INIT 
DO 60 I=1, NYEAR 
KYR = I 
NCT = NCT+1 
IF(NCT .EQ@. 2) GO TO 42 
CALL GGNORD(12, 13, Z) 
GO TO 44 
Z(1) = Z(2) 
NCT = 0 
TGI=TCII+SD#Z(1) 
TGI = TGII 
CALL SELECT 
CALL KILL (KYR) 
IF (IPERT. NE. 0. AND. MOD(AYR, INT). EG@.0) CALL DISTRB(KYR) 
IF (IPERT. NE. O. AND. KYR.EQ@. 1) CALL DISTRB(KYR) 
CALL BIRTH (KYR) 
CALL GROW(KYR) 
CALL OUTPUT (KYR, IPLOT) 
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(0120) 
(0121) 
(0122) 
(0129) 
(0124) 
(0125) 
(0126) 
(0127) 
(0128) 
(0129) 
(0130) 
(0131) 
(0132) 
(0133) 
(0134) 
(0135) 
(0136) 
(0137) 
(0138) 
(0139) 
(0140) 
(0141) 
(0142) 
(0143) 
(0144) 
(0145) 
(0146) 
(0147) 
(0148) 
(0149) 
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gaaNaaNAaAAaANa 


90 
60 
70 


80 


100 


110 


120 


150 


130 
140 


CALL PLTREE (KYR, IPLOT) 
CONTINUE 
CONT INUVE 
IF (IPLOT. NE. KALAST) GO TO 30 
ZNYR = KLAST 
DO 110 IV1=1,NYR1 
TIMEX = IV1i-1 
WRITE (7,80) TIMEX 
FORMAT (2E10. 3) 
DO 100 IV2=1,N3 
AC(TV1, IV2) = ACIV1, 1V2)/ZNYR 
WRITE (7,80) AC(IV1, IV2) 
CONTINUE 
CONT INUVE 


WRITE INFORMATION NEEDED FOR THE PLOT TO FILE TAPES. 
TAPES MAY BE EXAMINED OR THE PLOT PROGRAM APLOT MAY 
BE RUN. 


WRITE (8,120) NYEAR, NSPEC 
FORMAT (218) 
WRITE (8,150) ((SPNAM(CTI, J), 1=1,5), J=1, NSPEC) 
FORMAT (5A4) 
NNY = NYEAR+1 
DO 140 I=1, NNY 
WRITE (8,130) (ACI, J), J=1, NSPEC) 
FORMAT (15F8. 3) 
CONTINUE 
CALL EXIT 
END 


(0001 ) 
(0002) 
(0003) 
(0004) 
(0005) 
(0006) 
(0007) 
(0008) 
(0009) 
(0010) 
(0011) 
(0012) 
(0013) 
(0014) 
(0015) 
(0016) 
(0017) 
(0018) 
(0019) 
(0020) 
(0021) 
(0022) 
(0023) 
(0024) 
(0025) 
(0026) 
(0027) 
(0028) 
(0029) 
(0030) 
(0031) 
(0032) 
(0033) 
(0034) 
(0035) 
(0036) 
(0037 > 
(0038) 
(0039) 
(0040) 
(0041) 
(0042) 
(0043) 
(0044) 
(0045) 
(0046) 
(0047) 
(0048) 
(0049) 
(0050) 
(0051) 
(0052) 
(0053) 
(0054) 
(0055) 
(0056) 
(0057) 
(0058) 
(0059) 


cas 
cas 


gaananaannnannn 


eeokehehehenl| 


qgaangnaan 
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10 


20 


30 


SUBROUTINE ADREG(KYR) 
COMMON /FOREST/ NTREES(50), DBH( 700), IAGE (1000), KSPRT (50), NEWTR(50) 
1, SUMLA( 1000), NEW( 50), NCODE (700), ISEL (50) 


COMMON /PARAM/ AAA(100, 3), DMAX(50), DMIN(50), B39(50), B2(50), ITOL(50) 
COMMON /PARAM/ AAA(B8, 100), DMAX (50), DMIN( 50), B3(50), B2(50), ITOL(50) 
1, AGEMX (50), @(50), SPRTND( 50), SPRTMN(50), SPRTMX(50), KTIME(S0), Ai (50) 
2, A2(50), BB( 50), BC(50),C(50), WMIN( 50), IREG(S50, 4), DBHMX(50),B1(50), 
3B4(50), HMAX (S50), DRANGE (50) 
COMMON /CONST/ NSPEC, SOILG, PMS, TGI, TGII, SD, IGEOG 
COMMON /RAN/ YFL 
COMMON /DEAD/ NOGRO( 700) 
COMMON /COUNT/ NTOT, NYEAR 
COMMON /SEED/ USEED(13) 
INTEGER USEED 
REAL LA 
b Onder ADVANCE REGENERATION SUBROUTINE (ADREG) DIFFERS FROM GROW IN THAT 
aie: ole ONLY SHADE TOLERANT CONIFERS ARE ABLE TO GROW. THIS ROUTINE IS 
di Sido t CALLED AFTER WINDTHROWS 
bode EACH TREE IS REQUIRED TO ADD A 1.0 MM GROWTH RING EACH YEAR 
OTHERWISE THAT TREE IS SUBJECTED TO A POTENTIAL FOR SLOW- 
GROWTH RELATED MORTALITY (SUBROUTINE KILL) 


TINC = .1 


Bhat onic SUM LEAF AREA OF ALL TREES THAT ARE OF APPROXIMATELY 
Secale THE SAME HEIGHT 
DO 10 I=1, 700 
NOGRO(I) = O 
DO 20 I=1, 1000 
SUMLA(I) = O. 
SBIO = O. 
IF (NTOT .GT. 700) CALL ERR 
DO 50 K=1, NTOT 
I = NCODE(K) 
Suen: CALCULATE STAND BIOMASS, FOLIAGE BIOMASS IN KG, LEAF AREA 
INDEX 


LA IS THE LEAF AREA (8@. M) FOR PSME FROM SAPWOOD AREA 


LA=(256. #ALOG(DBH(K) )-639. )#2. 5 

IF (I. NE.6) GO TO 30 
HT=137. +B2(6) *#DBH(K)-B3(6) *DBH(K) ##2 
SBIO=SB10+. 0O2+2. O9*#DBH(K) ##2%(HT/100. )/100. -. 0015*#(DBH(K) #*#2% 

1 HT/100. )##2 

FBIO#3. 2+1. B9#ALOG (DBH(K) ) 

60 TO 31 
SBIO = SBIO+EXP(Ai(1)+BB(1)#ALOG(DBH(K) >) ) 

FBIO=EXP (A2(1)+BC(1)#ALOG(DBH(K) ) > 

IF (DBH(K). GE. 50. )FBIO=LA#(1/C (16) ) #DRANGE (I) 
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(0060) 35 SBIO = SBIO+FBIO 


(OOGT) Ca tes ear 

(OOG2)F7 1G err. HEIGHT PROFILE IS CALCULATED IN .1 METER UNITS 
(OOG3)P tC ee aaen 

(0064) 40 HIHT = (137+B2(1)*DBH(K)-B3(1)#DBH(K)##2)/10. +1. 
(0065) IF (HIHT. GT. 1000) GO TO 150 

(0066) IHT#HIHT 

(0067 ) IF(I.€@.5.0R.1.E@.6.0R. 1. EQ. 7.O0OR.1.E@.17.0R.1.E@.21)60 TO 45 
(0068) SUMLACIHT) = SUMLACIHT)+C (1) #FBIO/S5000. 

(0069) GO TO SO 

(0070) 45 SUMLA(IHT)=SUMLACIHT) + C(I) #FBIO/4000. 

(0071) SO CONTINUE 

(0072) DO 60 J=1,999 

(0073) Ji = 1000-J 

(0074) SUMLA(J1) = SUMLA(J1)+SUMLA(J1+1) 

(0075) 60 CONTINUE 

(OG7G)) C7 ene 

GOO77)| KC Seca CALCULATE AMOUNT OF GROWTH FOR EACH TREE 

(OO78D: SC” Shenae 

(0079) NG = 1 

(0080 ) DO 140 J#=1,NTOT 

(0081) I = NCODE(J) 

(0082) IF (ITOL(I) .NE. 5) GO TO 140 

COCES) Ne Co ere 

COOB4)) Gr Se raea CALCULATE MOISTURE —- GROWTH INDEX. 

(OCB8S)) Ce: 

(0086) GO TO (80, 90, 80. 80. 90, 70, 90, 90, 70, 90, BO, 90, BO, 90, 70, 90, 90, BO, 80 
(0087 ) 1 » 80,90), I 

COOSE) Re Ge lee MOIST SITE SPECIES 

(0089) 70 V = WMIN(I)/2. 

(0090) WATDEX = (PMS+2. )#(WMIN(I)—PMS) ##V/ (2. #WMIN( I) #*V) 
(0091) GO TO 100 

COOS2) TC yee INTERMEDIATE SITE SPECIES 

(0093) 80 Vi= (WMIN(T)-1. 079. 

(0094) WATDEX = (PMS+8. )#(WMIN(I)—-PMS)##V/ (9. #(WMIN(I)—-1. )##V) 
(0095) @0 TO 100 

GCOOS6) SC tesa DRY SITE SPECIES 

(0097 ) 90 Vi= (WMIN(I)—-5. 2/15. 

(0098) WATDEX = (PMS+10. )#(WMIN(I)—-PMS) ##V/(15. #(WMIN(I)—-5. )##V) 
(0099 ) 100 HT = 137+B2(1)#DBH(J)-B3(1)#DBH( J) ##2 

(0100) 102 CONT INUE 

(0101) IHT = HT/10. +2. 

(0102) 104 CONT INUE 

(0103) SLAR = SUMLACIHT) 

(0104) 4L = PHI#EXP(—-SLAR#. 26) 

(0105) LA=(256. *ALOG(DBH( J) )-639. )#*2.5 

(0106) FBIO@EXP (A2(1)+BC(1)#ALOG(DBH(U) ) >) 

(0107) IF (1. EG. 6)FBIO#3. 20+1. 8B9#ALOC (DBH( J) ) 

(0108) IF (DBH(J). GE. 50.) FBIO#@LA#¥(1. /C (16) )#*#DRANGE (I) 
(0109) 120 A = (1. -(DBH(J)#HT) / (DBHMX (I) *HMAX (I) )) 

(0110) B= DBH(J)#B1(1)#FBIO##B4(1>) 

(0111) CC = (274+3. #B2(1)#DBH( J) ##2-4. #B3(1)#DBH( J) ##*3) 
(0112) D= (1. -SBIO/SOILQ) #4 

(0113) TEGD = (TCI-DMIN(I))#(DMAX(1)-TGI)/(DMAX(1)-DMIN(I) )##2 
(0114) DNC = (A * B / CC) * D * TEGD #WATDEX 

(0115) DINC = (1. -EXP(-32. 846% (AL-O. 025) ) )*DNC 

(0116) IF (DBH(J) .LT. DBHMX(I)) GO TO 121 

(0117) DINC = 0.0 

(0118) 121 CONT INUVE 
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qgaaaanNn 


abet CHECK INCREMENT LESS THAN 1.0 MM REQUIRED GROWTH 


peendiaps WRITE OUT DIAMETER INCREMENT FOR PSME 


IF (I .NE. 16) 6O TO 125 
IF (DBH(J) .LE. 86. .OR. DBH(J) .@T. 165.) 60 TO 125 
WRITE (9,122) DINC, KYR 
122 FORMAT (1X, F9. 4, 3X. 15) 
125 CONTINUE 
IF (DINC. GE. TINC) GO TO 130 
IF (DINC.LT.O.) DINC = O. 
NOGRO(NG) = J 
NG = NG+1 
130 DBH(J) = DBH(J)+DINC 
140 CONTINUE 
RETURN 
150 WRITE (6, 160) 
160 FORMAT (19H1 IHT EXCEEDED 1000) 
STOP 
END 
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SUBROUTINE BIRTH (KYR) 

COMMON /FOREST/ NTREES(S0), DBH(700), IAGE( 1000), KSPRT (50), NEWTR( SO) 
1, SUMLA(1000), NEW(S0), NCODE (700), ISEL(50) 

COMMON /PARAM/ AAA(8, 100), DMAX(50), DMIN( 50), B3(50), B2(50), ITOL( 50) 
1, AGEMX (50), 6(50), SPRTND( 50), SPRTMN( 50), SPRTMX (50), ATIME(50), A1(50) 
2, A2(50), BB( 50), BC(50),C(50), WMIN( 50), IREG(S50, 4), DBHMX(50),B1(50), 
3B4 (50), HMAX (50), DRANGE (50) 

COMMON /CONST/ NSPEC, SOILG, PMS. TGI, TGII, SD, IGEOG 

COMMON /RAN/ YFL 

COMMON /DEAD/ NOGRO( 700) 

COMMON /COUNT/ NTOT. NYEAR, IPERT, INT, IN, PMORT 

COMMON /SEED/ USEED(13) 

INTEGER USEED 

REAL LA 

Bn Don SAPLINGS ENTER THE PLOT AT AVERAGE SIZE OF 10.0 CM DBH 

IF(KYR.LE. 20) GO TO 3 
IF (NTOT. GT. 200. . AND. MOD(KYR, 10) .E@. 0) GO TO 120 
CONT INUE 

SIZE = 10.0 

Srey SELECT A SPECIES FROM 1 TO 3 TIMES 

....ADD UP NUMBER OF ABAM AND TSHE OF SIZE 10 TO 20 CM SO THAT 

A RESTRICION CAN BE PUT ON THE NUMBER OF BIRTHS. 
NABAM = 0 
NTSHE = 0 
DO 9 K = 1, NTOT 
KK=NCODE (KK) 
IF (KK.NE.1) GO TO 5 
IF (DBH(K). GE. 20. .OR. DBH(K) .LT. 10.) 6CO TO 5 
NABAM = NABAM + 1 
CONTINUE 
IF (KK.NE. 19) GO TO 8 
IF (DBH(K). GE. 20. .OR. DBH(K) .LT. 10.) GO TO 8 
NTSHE = NTSHE + 1 
CONT INUVE 
CONT INUE 
NPLANT = 3. #RANF (6) +1. 


DO 100 JK=1, NPLANT 
FBIO = O. 
SED Ac CALCULATE BIOMASS (LEAF MASS) FOR EACH SPECIES 
IF (NTOT.E@.0) GO TO 40 
DO 30 K=1, NTOT 
KA = NCODE(K) 
LA=(256. O#ALODG@ (DBH(K) )-639. ) 
IF (AK. E@. 5. OR. KK. EQ. 7. OR. KK. EG. 17. OR. KK. EG. 21) CO TO 25 
FBIO=@EXP (A2(KK)+BC (KK) #ALOG(DBH(K))) #C(KK)/2.5 +FBIO 
IF (DBH(K). OE. 50.) FBIO=FBIO+LA#DRANGE (KK) /2. 5 
60 TO 30 


FBIO=FBIO +EXP(A2(KK)+BC (KK) #ALOG(DBH(K))) *#C(KK)/2. 0 
IF (DBH(K). GE. SO.) FBIO=FBIO+LA#DRANOGE (KK) /2. 0 

CONT INUVE 

NW = 0 


XLAI = FBIO/2000. 


= | 


(0060) 
(0061) 
(0062) 
(0063) 
(0064) 
(0065) 
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(0068) 
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(0084) 
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(0096) 
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(0098) 
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(0107) 
(0108) 
(0109) 
(0110) 
(0111) 
(0112) 
(0113) 
(0114) 
(0115) 
(0116) 
(0117) 
(0118) 


(0119) 
(0120) 
(0121) 
(0122) 
(0123) 
(0124) 
(0125) 


aan anNn gnaNaaaNNn aan 


aan 


aan 


DO 70 J=1, NSPEC 
IF (J. EQ@.1 .AND. NABAM .@GT. 200) @O TO 70 
IF (J.E@.19 . AND. NTSHE .GT. 100) 60 TO 70 
IF (ISEL(J).NE.1) CO TO 70 
IF (XLAI.LT. 3.) 6O TO 60 
@O TO (50, 50, 50, 50, 50, 70, 50, 50, 50, 70, 50, 50, 50, 70, 50, 50, 70 
» 30, 30, SO, 30), J 
IF (XLAI.LT. 6.) GO TO 60 
GO TO (54, 54, 34, 70, 54, 54, 54, 54, 54, 54, 54, 54, 54, 54, 54, 70, 
54,54, 54,54,70), J 
IF (XLAI .LT. 8.) GO TO 60 
@O TO (58, 70. 70, 58. 70, 58. 70. 58, 58, 58, 70, 70, 70, 58, 70, 58, 
58, 58, 58,58, 58). J 
IF (XLAI .LT. 10.) 60 TO 60 
60 TO (60, &O, 40, 60, 60, 60, 60, 70, 70, 60, 60, 40, 60, 60, SO, 40, 40, 
70,60, 70,60), J 
NW = NW+1 
NEWTR(NW) = J 
CONT INVE 


NW = NW*RANF(10)+1.0 
NSP = NEWTR (NW) 


__.. PLANT RANDOM NUMBER OF SEEDLINGS 


DO 80 J=1, MPLANT 
NTOT = NTOT+1 
IF (NTOT. GT. 700) CALL ERR 
IAGE(NTOT) = 0 
NCODE(NTOT) = NSP 


|... CALCULATE DBH FOR SEEDLINGS 


DBH(NTOT) = SIZE+5. #(1. O-RANF (11) )##3 
NTREES(NSP) = NTREES(NSP)+1 

CONT INVE 

CONT INVE 

CALL SPROUT 


_... BIRTH CONTINUED UNTIL LAI .@E. 1. 


XLAI = FBIO/2000. 
IF (XLAI.LT.1.)> 60 TO 10 


100 CONTINUE 


.... INCREMENT AGES 


DO 110 I=1, NTOT 


IAGE(I) = IAGE(I)+1 


110 CONTINUE 
120 CONTINUE 


RETURN 
END 


43 


(0001) 
(0002) 
(0003) 
(0004) 
(0005) 
(0006) 
(0007 ) 
(0008) 
(0009) 
(0010) 
(0011) 
(0012) 
(0013) 
(0014) 
(0015) 
(0016) 
(0017) 
(0018) 
(0019) 
(0020) 
(0021 ) 
(0022) 
(0023) 
(0024) 
(0025) 
(0025) 
(0027 ) 
(0028) 
(0029) 
(0030) 
(0031) 
(0032) 
(0033) 
(0034) 
(0035) 
(0036) 
(0037) 
(0038) 
(0039) 
(0040) 
(0041) 
(0042) 
(0043) 
(0044) 
(0045) 
(0046) 
(0047) 
(0048) 
(0049) 
(0050) 
(0051) 
(0052) 
(0053) 
(0054) 
(0055) 
(0056) 
(0057) 
(0058) 
(0059) 


44 


c@s 


ieiigity) 


gaan 


30 
40 


45 
50 


60 
65 


68 


OO OUDICNA AO OC. arg 


SUBROUTINE DISTRB(KYR) 
COMMON /FOREST/ NTREES(50), DBH( 700), IAGE(1000), KSPRT(50), 
1 NEWTR( 50), SUMLA( 1000), NEW(S50), NCODE(700), ISEL(50) 


COMMON /PARAM/ AAA(8, 100), DMAX (50), DMIN( 50), B3(50), B2(50), ITOL(S0) 
1, AGEMX (50), (50), SPRTND( 50), SPRTMN( 50), SPRTMX (50), ATIME (50), A1(50) 
2,A2(50), BB( 50), BC(50),C(50), WMIN( S50), IREG( 50, 4), DBHMX(50), B1(50), 


3B4(50)., HMAX (50), DRANGE (50) 


COMMON /CONST/ NSPEC,. SOILG, PMS, TGI, TGII., SD, IGEOG 
COMMON /COUNT/ NTOT, NYEAR, IPERT, INT, IN, PMORT, FI 
COMMON /SEED/ USEED(13) 

INTEGER USEED 

IF (IPERT. EQ@. 0) RETURN 

IF (IPERT.EQ@. 1) WRITE(6, 16) INT 

IF (IPERT. EQ. 2) WRITE(6, 17) INT 

IF (IPERT. EG. 3) WRITE(6, 18) INT 

FORMAT (1HO, 11H FIRE EVERY, 15, 46H YEARS) 

FORMAT (1HO, 17H WINDSTORMS EVERY, 15, 6H YEARS) 
FORMAT (1HO, 13H CLEARCUTS ON, I5, 9H ROTATION) 
IF(IPERT .NE. 1) GO TO 65 


FIRE 
DO 60 K=1, NTOT 


I=NCODE (K) 
P=0.0 


IF (I. NE. 17.0R. 1. NE. 12. OR. I. NE. 14. OR. I. NE. 16. OR. I. NE.8)GO TO 30 


FIRE TOLERANT SPECIES HAVE THE FOLLOWING PROBABILITY 
OF DYING IN A FIRE 


P=(EXP(-. 0O0255#*DBH(K) ) ) *FI 
GO TO 40 

P=(EXP(—. OOOS3#DBH(K) ) )#FI 
IF (RANF(2) .GT. P) GO TO 60 
NTREES(1I)=NTREES(I)-1 


IF (DBH(K). GT. SPRTMN( I). AND. DBH(K). LT. SPRTMX( 1) )ASPRT(1I)=—-1 


IF (PMORT .NE. 9) GO TO 50 
WRITE (9,45) NCODE(K), DBH(K),KYR 
FORMAT (3X,13,F6. 2, 14) 
CONT INVE 
DBH(K)=-1 
CONT INUE 
CONT INVE 
IF (IPERT .NE. 2) 6O TO 98 


LIGHT WINDSTORMS BASED ON THE TAR CREEK DATA (OF J. FRANKLIN) 


DO 80 K =1,NTOT 


I=NCODE (K) 
LFA CISEGI1T6 ORS IREGA 13) SOR. sb Qst2 
1 .OR. 1.E@.14 .OR. I.EQ@. 4) 60 TO 70 


TE CIS EQS) SOR REG 6) JOR] hI EGY 7a OR ale Gilet OR: 


1 I1.—€@.18 .OR. I.E@.21 .OR. I.EQ@.9 .OR. I.EQ@.8) 60 TO 74 


IF (DBH(K) .@T. 60.) GO TO 68 
P = .15 -.001 * DBH(K) 

GO TO 79 
Pae=s 107, 

GO TO 79 


(0060) 
(0061) 
(0062) 
(0063) 
(0064) 
(0065) 
(0066) 
(0067 ) 
(0068) 
(0069) 
(0070) 
(0071) 
(0072) 
(0073) 
(0074) 
(0075) 
(0076) 
(0077) 
(0078) 
(0079) 
(0080 ) 
(0081 } 
(0082) 
(0083) 
(0084) 
(0085) 
(0086) 
(0087 ) 
(0088) 
(0089) 
(0090) 
(0091) 
(0092) 
(0093) 
(0094) 
(0095) 
(0096) 
(0097 ) 
(0098) 
(0099) 
(0100) 
(0101) 
(0102) 
(0103) 
(0104) 
(0105) 
(0106) 
(0107) 
(0108) 
(0109) 
(0110) 
(O111) 
(0112) 
(0113) 
(0114) 
(0115) 
(0116) 
(0117) 
(0118) 


G. 70) IF (DBH(K) .GT. 100.) GO TO 72 
Cc P = .13 - .001 * DBH(K) 
C GO TO 79 
Cc 72 P= .03 
Cc GO TO 79 
Cc 74 IF (DBH(K) .GT. 100.) GO TO 78 
C P = .30 -.02 # DBH(K) 
Cc co TO 79 
Cc 78 P= .18 
Cc 
C79, IF(RANF(3) .@T.P) 60 TO 80 
C NTREES(I) = NTREES(I)-1 
Cc IF (DBH(K). GT. SPRTMN(I). AND. DBH(K). LT. SPRTMX(I)) 
C 1 ASPRT(1I)=-1 
C DBH(K)=-1 
80 CONTINUE 


90 CONTINUE 
CATASTROPHIC WIND 
DO 96 K=1,NTOT 
Cc SPRUCE & HEMLOCKS ARE WIND INTOLERANT 
I = NCODE(K) 
IF (1.E€@.11 .OR. I .E@.19 .OR. I1.EQ. 20) 
WRITE(6,91) DBH(K),K 


91 FORMAT ¢ F10. 3,15) 
P=(1-(EXP(-. O0255#DBH(K) ))+0. 6) #FI 
60 TO 94 
92 P=(1-(EXP(-. OOOS3*#DBH(K) ))+0. 9) #FI 
94 IF (RANF(3).GT.P) GO TO 96 


60 TO 92 


IF (DBH(K). @T. SPRTMN(I) . AND. DBH(K).LT. SPRTMX(I)) 


1 KSPRT(1I)=-1 
IF (PMORT .NE. 9) GO TO 95 
WRITE (9,45) NCODE(K), DBH(K), KYR 
95 CONT INUE 
DBH(K)=-1 
96 CONT INUE 
98 IF (IPERT. NE. 3) GO TO 110 


CLEARCUTS 


aan 


DO 100 K=1, NTOT 
NTREES (I) =NTREES(I)-1 
IF (PMORT .NE. 9) GO TO 99 
WRITE (9,45) NCOQDE(K), DBH(K), KYR 
99 CONT INVE 
DBH(K)=-1 
100 CONT INUVE 
110 CONTINUE 
148 K1=0 
DO 150 K=1, NTOT 
IF(DBH(K).LT.O.) GO TO 150 
KA1izKi+li 
DBH(K1)=DBH(K) 
IAGE(K1)=IAGE(K) 
NCODE (K1 ) =NCODE (K) 
150 CONTINUE 
160 CONT INUE 
NTOT=K1 
IF (IPERT .NE. 3) GO TO 160 
Cc AFTER CLEARCUT PLANT PSEUDOTSUGA MENZIESII 


(OTHER SPECIES COULD BE 
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(0124) 
(0125) 
(0126) 
(0127) 
(0128) 
(0129) 
(0130) 
(0131) 
(0132) 
(0133) 
(0134) 
(0135) 
(0136) 
(0137) 
(0138) 
(0139) 
(0140) 
(0141) 
(0142) 
(0143) 
(0144) 
(0145) 
(0146) 
(0147) 
(0148) 
(0149) 
(0150) 
(0151) 


(0001) 
(0002) 
(0003) 
(0004) 
(0005) 


(0001) 
(0002) 
(0003) 
(0004) 
(0005) 
(0006) 
(0007 ) 
(0008) 
(0009) 
(0010) 
(0011) 
(0012) 
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C 
Cc 


170 


220 


10 


PLANTED BY MAKING CHANGES IN THE "170" LOOP) 


DO 170 J = 1,50 

NTOT=NTOT + 1 

TAGE (NTOT)=0 

NCODE (NTOT) #16 

DBH(NTOT)= 10. + 5. #(1. O-RANF (11) )##3 

NTREES(16) = NTREES(16) + 1 
CONTINUE 
INCREASE YEAR COUNT BY 30 TO MAKE UP LOSS SINCE THESE TREES 

WERE INTRODUCED AT 10 CM DBH AFTER THE CLEARCUT. 


KYR=KYR+30 


CAUSE ADVANCE REGENERATION BY CALLING SUBROUTINE ADREG AFTER WIND 
IF (IPERT. NE. 2) GO TO 200 


DO 190 IT=1,3 
CALL ADREG(KYR) 
CONT INVE 
CONTINUE 
SINCE TREES ARE INTRODUCED AT 10 CM DBH, LET MODEL GROW TREES 
FOR 30 YEARS TO MAKE UP LAG IN GROWTH. 


IF (IPERT.NE.1 .OR. IPERT .NE. 2) GO TO 220 
DO 210 I=1,30 
CALL KILL(KYR?) 
CALL GROW(KYR) 
KYR=KYR+1 
CONTINUE 
RETURN 
END 


SUBROUTINE ERR 

WRITE (6,10) 

FORMAT (37H1 THE NUMBER OF TREES HAS EXCEEDED 800) 
STOP 

END 


SUBROUTINE GGNORD (NSEED1, NSEED2, Z) 
DIMENSION Z(1) 

DATA PI2/0. 62831853E01/ 

K = 0 

Al = RANF (NSEED1) 

A2 = RANF (NSEED2) 

K = K+l1 

Z(K) = SQRT(-. 2EO1#ALOC(A1) )*#SIN(PI2#A2) 
K = K+1 

Z(K) = SQRT(-0. 2E01#ALOG(A1) )#COS(PI2*A2) 
RETURN 

END 
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SUBROUTINE GROW(KYR ) 
COMMON /FOREST/ NTREES(50), DBH(700), IAGE( 1000), KSPRT (50), NEWTR(50) 
1, SUMLA( 1000), NEW( 50), NCODE( 700), ISEL(S50) 


COMMON /PARAM/ AAA( 100, 3), DMAX(50), DMIN( 50), B3(50), B2(50), ITOL(50) 
COMMON /PARAM/ AAA(B, 100), DMAX (50), DMIN( 50), B3(50), B2(50), ITOL (50) 
1, AGEMX (50), 6(50), SPRTND( 50), SPRTMN( 50), SPRTMX (50), ATIME(50), A1 (50) 
2, A2(50), BB( 50), BC(50),C(50), WMIN( SO), IREG(50, 4), DBHMX (50), B1(50), 
3B4 (30), HMAX (50), DRANGE (50) 
COMMON /CONST/ NSPEC, SOILQG, PMS, TGI, TGII, SD, IGEOG 
COMMON /RAN/ YFL 
COMMON /DEAD/ NOGRO(700) 
COMMON /COUNT/ NTOT, NYEAR 
COMMON /SEED/ USEED(13) 
INTEGER USEED 
REAL LA 
Marie EACH TREE IS REQUIRED TO ADD A 1.0 MM GROWTH RING EACH YEAR 
OTHERWISE THAT TREE IS SUBJECTED TO A POTENTIAL FOR SLOW- 
GROWTH RELATED MORTALITY (SUBROUTINE KILL) 


TINC = .1 


ete SUM LEAF AREA OF ALL TREES THAT ARE OF APPROXIMATELY 
Ame ee THE SAME HEIGHT 
DO 10 I=1, 700 
NOGRO(I) = 0 
DO 20 I=1, 1000 
SUMLA(I) = O. 
SBIO = O. 
IF (NTOT .GT. 700) CALL ERR 
BO 50 K=1, NTOT 
I = NCODE(K) 
Son ek CALCULATE STAND BIOMASS, FOLIAGE BIOMASS IN KG, LEAF AREA 
INDEX 
LA IS THE LEAF AREA (SG. M) FOR PSME FROM SAPWOOD AREA 
LA=(256. #ALOG(DBH(K) )-639. )#2. 5 
IF (1I.NE.6) 60 TO 30 
HT=137. +B2(6) #*DBH(K)—-B3( 46) #DBH(K) ##2 
SBIO=SBIO+. 0O2+2. O9#DBH(K) #¥#2*%(HT/100. )/100. —. 0015# (DBH(K) ##2% 
1 HT/100. )*##2 
FBIO#=3. 2+1. B9#ALOG(DBH(K) ) 
e0 TO 31 
SBIO = SBIO+EXP(Ai1(1)+BB(1I)#ALOCG(DBH(K) ) ) 
FBIO=EXP (A2(1)+BC(1)#ALOG(DBH(K) ) ) 
IF (DBH(K). GE. 50. )FBIO=LA#(1/C (16) ) *DRANGE (1) 
SBIO = SBIO+FBIO 


Ser HEIGHT PROFILE IS CALCULATED IN .1 METER UNITS 


HIHT = (137+B2(1)*DBH(K)-B3(1)#DBH(K) *#2)/10. +1. 
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(0060) IF (HIHT. GT. 1000) GO TO 150 


(0061) IHT#HIHT 

(0062) IF(I.€Q@.5.0OR.1.E€Q@. 6 OR. 1.EQ@. 7.OR.1.EQ@.17.0R. 1. EQ. 21)60 TO 45 
(0063) SUMLACIHT) = SUMLA(IHT)+C(1)#FBIO/5000. 

(0064) eO TO 50 

(0065) 45 BUMLA(CIHT)=SUMLACIHT) + C(I) #FBIO/4000. 

(0066) 5O CONTINUE 

(0067) DO 60 J=1,999 

(0068) J1 = 1000-J 

(0069 ) SUMLA(J1) = SUMLA(J1)+SUMLA(J1+1) 

(0070) 60 CONTINUE 

COOZAS) Ga halwatons 

€OO72)) Cee CALCULATE AMOUNT OF GROWTH FOR EACH TREE 

COO73)) iC en 

(0074) NG = 1 

(0075) DO 140 J=1,NTOT 

(0076) I = NCODE(J) 

(0077) {Ce ee 

(OO7B)— Cy ceicine. CALCULATE MOISTURE —- GROWTH INDEX. 

COO7D MGA it Wacneror 

(0080 ) GO TO (80. 90, BO. BO, 90, 70, 90, 90. 70, 90, BO, 90, BO, 90, 70, 90, 90, 80, 80 
(0081) 1 »80,90), I 

(OOB82) 5 Ca cree MOIST SITE SPECIES 

(0083) 70 V = WMIN(I) 7/2. 

(0084 ) WRITE (6,72) I, WMIN(I),Y 

(00e5) 72 FORMAT (6H I= ,13,6H WMIN= ,F6.2,6H V= ,F6. 2) 
(0086) WATDEX = (PMS+2. )#(WMIN(I)—-PMS)##V/ (2. #WMIN( I) ##V) 
(0087 ) WRITE (6,73) WATDEX 

(0088) 73 FORMAT ( 8H WATDEX= ,F6. 2) 

(0089) 60 TO 100 

COO9O)) Cy sere INTERMEDIATE SITE SPECIES 

(0091) BO V = (WMIN(I)-1. )/9. 

(0092) WATDEX = (PMS+8. )#(WMIN(I)—PMS)##V/(9. #(WMIN(I)—-1. )##V) 
(0093) 60 TO 100 

COOP 4) EC ee ein DRY SITE SPECIES 

(0095) 90 Vi = (WMIN(I)-5. 2/15 

(0096) WATDEX = (PMS+10. )#(WMIN(1)-PMS)##V/(15. #(WMIN(I)—-5. )*#V) 
(0097) 100 HT = 137+B2(1)*DBH(J)-B3(1)#DBH( J) ##2 

(0098) 102 CONT INUE 

(0099) IHT = HT/10. +2. 

(0100) 104 CONT INUVE 

(0101) SLAR = SUMLA(IHT) 

(0102) AL = PHI#EXP(-SLAR*. 26) 

(0103) LA=(256. #ALOG(DBH(J) )-639. )#2.5 

(0104) FBIO@EXP (A2(1)+BC(1I)#ALOG(DBH( J) )) 

(0105) IF(I. EQ. 6)FBIO=3. 20+1. 8B9*ALOG(DBH( JU) ) 

(0106) IF (DBH(J). GE. 50. ) FBIO=LA*(1. /C (16) )#DRANGE (I> 
(0107) 120 A = (1.-(DBH(JU)#HT) / (DBHMX (1) #HMAX(1))) 

(0108) B= DBH(J)#B1(1)#FBIO##B4(1) 

(0109) CC = (274+3. #B2(1)#DBH( JU) ##2-4. #B3(1)*DBH( JU) ##3) 
(0110) D= (1.-SBIO/SOILQ@) #4 

(0111) TEGD = (TGI-DMIN(I) )#(DMAX(1I)-TGI)/(DMAX(1)-DMIN(I) ) *##*2 
(0112) DNC = (A * B / CC) * D * TEGD *#WATDEX 

(0113) DINC = (1. -EXP(-20. 467#(AL-O. 025) ) ) *DNC 

COTS EC ere 

GOULD SIG. nora 

COUUG Ca ee crear: CHECK FOR PS CURVE TYPE —- PONDEROSA PINE 

CO6197 2) as C ies, oie 


(0118) IF (ITOL(I).EQ@.1) DINC = (1. O-EXP(-23. 026*(AL-O. 30) ))#*DNC 
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(0119) IF(ITOL(I).EQ@.1 . AND. AL.LT.O.3) DINC=O. 


(OLNZO) (Gea 

(OTe) Ce ae CHECK FOR PS CURVE TYPE — HARDWOODS 

CON22)) Ce Beet 

(0123) IF (ITOL(I).EQ@.2) DINC = (1. O-EXP(-9. 695#(AL-\). 025) ))#DNC 
CONZE NSS GY sy 

COMA SI Cis Tele arene CHECK FOR PS CURVE TYPE — TOLERANT CONIFERS 
(0126) IF (ITOL(I). €Q@.4) DINC = (1. -EXP(-32. 8B46%(AL-O. 025) ) )*#DNC 
(0127) IF (ITOL(I).EQ@.5) DINC = (1. -EXP(-40. O#(AL-O. 025) ) )#DNC 
(0128) IF (DBH(J) .LT. DBHMX(I)) GO TO 121 

(0129) DINC = 0.0 

(0130) 121 CONT INUE 

COUSIDT Ce err 

COLS209 1G ae CHECK INCREMENT LESS THAN 1.0 MM REQUIRED GROWTH 
COUSS)I GC ee 

(COVS4) (Co ae ae WRITE OUT DIAMETER INCREMENT FOR PSME 

(0135) C 

(0136) IF (I .NE. 16) 60 TO 125 

(0137) IF (DBH(J) .LE. 86. .OR. DBH(J) .GT. 165.) @0 TO 125 
(0138) IF (PMORT. NE.2) GO TO 125 

(0139) WRITE (9%,122) DINC, KYR 

(0140) 122 FORMAT (1X,F9. 4, 5X, 15) 

(0141) 125 CONTINUE 

(0142) IF (DINC. GE. TINC) GO TO 130 

(0143) IF (DINC.LT.0O.) DINC = O. 

(0144) NOGRO(NG) = J 

(0145) NG = NG+1 

(0146) 130 DBH(J) = DBH(J)+DINC 

(0147) 140 CONTINUE 

(0148) RETURN 

(0149) 150 WRITE (6, 160) 

(0150) 160 FORMAT (19H1 IHT EXCEEDED 1000) 

(0151) STOP 

(0152) END 

(0001 ) SUBROUTINE INIT 

(0002) COMMON /FOREST/ NTREES(50), DBH( 700), IAGE (1000), KSPRT (50), NEWTR(SO) 
(0003) 1, SUMLA( 1000), NEW(50), NCODE (700), ISEL (50) 

(0004) COMMON /CONST/ NSPEC, SOILG, PMS, TGI, TCII, SD, IGEOG 
(0005) COMMON /DEAD/ NOGRO( 700) 

(0006) DO 10 I=1, N65PEC 

(0007) NOGRO(I) = O 

(0008) ASPRT(I) = 1 

(0009) 10 CONTINUE 

(0010) NSPE1 = NSPEC+1 

(0011) DO 20 I=NSPE1, 700 

(0012) NOGRO(I) = O 

(0013) 20 CONTINUE 

(0014) Cas RND IS THE RANDOM NUMBER GENERATOR ON THE PRIME 

(0015) CALL R=RND(1) 

(0016) RETURN 


(0017) END 


(0001) SUBROUTINE INPUT 


(0002) COMMON /PARAM/ AAA(8, 100), DMAX(50), DMIN( 50), B3(50), B2(50), ITOL(50) 
(0003) 1, AGEMX (50), @(50), SPRTND(50), SPRTMN( 50), SPRTMX (50), KTIME( 50), A1(50) 
(0004) 2,A2(50), BB( 50), BC(50),C( 50), WMIN( SO), IREG( 50, 4), DBHMX(50), B1(50), 
(0005) 3B4(50), HMAX (50), DRANGE (50) 

(0006) COMMON /CONST/ NSPEC, SOILG, PMS, TGI, TGII, SD, IGEOG 

(0007 ) COMMON /COUNT/ NTOT, NYEAR, IPERT, INT. IN, PMORT, FI 

(0008 ) DIMENSION NSELCT(50) 

(0009) READ (5.10) NSPEC, (NSELCT(I), I=1, NSPEC) 

(0010) 10 FORMAT (4012) 

(0011) READ (5,20) TGII.,SD, PMS, IGEOG., NYEAR, IPERT, INT, PMORT. IN 

(0012) 20 FORMAT (FS. 1, 2F3.0, 11, 514) 

(0013) ¢ 

(OO Caen  e 
(0015) C NYEAR#=LENGTH OF THE MODEL RUN IN YEARS 

(0016) C IPERT= DISTURBANCE TYPE 

(0017) C¢ O NONE 

(0018) C 1 FIRE 

(0019) C 2 WINDTHROW 

(0020) C¢ 3 CLEARCUT 

(0021) C INT=INTERVAL OF DISTURBANCE IN YEARS 

(0022) C PMORT - IF EQUAL TO 

(0023) C 2 = PRINT OUT DIAMETER INCREMENT OF DOUGLAS FIR ON YEARLY BASIS 
(0024) C 7 = PRINT OUT STAND STRUCTURE AT YEAR 700 FOR INPUT INTO NEXT 
(0025) C RUN OF THE MODEL 

(0026) C 9 = PRINT OUT MORTALITY INFORMATION: SPECIES, DBH AND YEAR DIED 
(0027) C FOR INPUT INTO DECAY MODEL - SEE GRAHAM, R. 1982. PH.D. THESIS 
(0028) C AT OREGON STATE UNIVERSITY, P. 112. 

(0029) C 

(0030) C IN=READ IN STAND STRUCTURE IF NOT O 

(OO31) C = NUMBER OF TREES TO BE READ IN (THE DATA IS AT THE END OF THE 
(0032) C FILE SUCDAT> 

(0033) C 

SCOX ORE PPB). SSI SSO SSS 
(0035) READ (5,25) FI 

(0036) 25 FORMAT (F5. 2) 

(0037) C FI= FIRE OR WIND INTENSITY (O<rI<=1) 

(0038) C 

(0039) DO 40 J=1,NSPEC 

(0040) READ (5,30) (AAACI, J), 11,5), DMAX( J), DMIN( J), B3(J), B2( J), ITOL( J 
(0041) 1 }, AGEMX (J), SPRTND( J), SPRTMN( J), SPRTMX (J), KTIME (JU), NUM, A1( J), A2 
(0042) 2 (J), BB(J), BC( J), C( J), WMIN( J), (IREG(J, 1), I=1, 4), DBHMX(J), BI (J), 
(0043) 3 HMAX (J), B4( J), DRANGE (JU) 

(0044) 30 FORMAT (SA4 ,6xX,F4.0,1X, 2F4. 0,F5.0,11,F4. 0, 7X, F2. 0, 2F4. 0, 14, 3X, 12/ 
(0045) 1F7.4,F6.3,F7.4,F6. 4,1X,F4.1,F3.0,411/F4. 0,F5. 0,F5.0,F3. 1,F6. 3) 
(0046) 

(0047) 40 CONTINUE 

(0048) C 

(0049) C FOR “NONSTRESSED" SITES USE HEIGHT PARAMETERS DERIVED FROM 
(0050) C REGRESSION ANALYSIS. 

(0051) C 

(0052) IF(TGII .@GT. 60.0 . AND. PMS .LT. 15.) GO TO 45 

(0053) B2(16)=72. 34 

(0054) B3(16)=0. 22 

(0055) B2(19)275. 97 

(0056) B3(19)=0. 27 

(0057) 45 CONT INUVE 

(0058) WRITE (6, 50) 

(0059) SO FORMAT (1Hi, 1H . 6X, 7HSPECIES, 11X, 11HDEGREE DAYS, 2X, 11HHEIGHT PARA, 
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(0060) 
(0061) 
(0062) 
(0063) 
(0064) 
(0065) 
(0066) 
(0067 ) 
(0068) 
(0069) 
(0070) 
(0071) 
(0072) 
(0073) 
(0074) 
(0075) 
(0076) 
(0077) 
(0078) 
(0079) 
(0080) 
(0081) 
(0082) 
(0083) 
(0084) 
(0085) 


12X, 4HITOL, 2X, SHAGEMX, 2X, GHSPRTND, /, 26X, GHMAX, 3X, GHMIN, 4X, 2HBS, 4X, 
22HB2/ ) 
DO 70 I=1, NSPEC 
WRITE (6,60) I, (AAACJ, 1), J=1,5), DMAX(I), DMIN( IT), B3(1), B2(T), 
1 ITOL (1), AGEMX (1), SPRTND(T) 
60 FORMAT (1H ,12,1X, 5A4 , 2X, F4. 0, 3X, F4. 0, 1X, F5. 2,1X,F6. 2, 3X, 12, 3X, 
1FS. 0. 4X, F3. 0) 
70 CONTINUE 
WRITE (6, 80) 
BO FORMAT (////, 1H . 6X, 7HSPECIES, 11X, 6HSPRTMN, 2X, GHSPRTMX, 2X, SHATIME, 
16X, QHA1, 7X, 2HA2, 7X, 2HBB, 7X, 2HBC, 7X, 1HC, 6X, 4HWMIN/ ) 
DO 100 I=1, NSPEC 
WRITE (6,90) I, (AAACJ, I), J=1,5), SPRTMN(I), SPRTMX(1),KTIME(TI), Al 
1 (I), A2(1), BBCI), BCCI), C(I), WMIN(T) 
90 FORMAT (1H ,12,1X, 5A4 .2X,F4. 0, 4X, F4. 0, 2X, 16, 6(2X, F7. 4)) 
100 CONTINUE 
WRITE (6,110) 
110 FORMAT (1H1.,1H . 6X, 7HSPECIES, 16X, 4HIREG, 6X, SHDBHMX, 6X, 2HB1, 6X, 
14HHMAX, 8X, 2HB4, BX, GHDRANGE, /, 28X, 1H1, 2X. 1H2, 2X, 1H3, 2X, 1H4/) 
DO 130 I=1, NSPEC 
WRITE (6,120) I, (AAACJ, I). Jel, 5), (IREG(I, J), J=1, 4), DBHMX(1), Bi 
1 (1), HMAX(1), B4( 1). DRANGE(T) 
120 FORMAT (1H ,I2, 1X, 5A4 , 2X, 413, 3X, F5. 0. 3(2X, FB. 1), 2X, Fe. 3) 
136 CONTINUE 
RETURN 
END 
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SUBROUTINE KILL (KYR) 

COMMON /FOREST/ NTREES( 50), DBH(700), IAGE( 1000), KSPRT( 50), NEWTR(S0) 
1, SUMLA( 1000), NEW(50), NCODE(700), ISEL (50) 

COMMON /PARAM/ AAA(8, 100), DMAX(50), DMIN( 50), B3(50), B2(50), ITOL( SO) 
1, AGEMX (50), @(50), SPRTND(50), SPRTMN( 50), SPRTMX(50), KTIME(50), A1 (50) 
2, A2(50), BB( 50), BC(50),C(50), WMIN( 50), IREG(50, 4), DBHMX(50), Bi(50), 


3B4(50), HMAX (50), DRANGE (50) 


COMMON /CONST/ NSPEC, SOILG, PMS, TGI, TGII, SD, IGEOC 
COMMON /DEAD/ NOGRO( 700) 

COMMON /COUNT/ NTOT. NYEAR, IPERT, INT. IN, PMORT 
COMMON /SEED/ USEED(13) 

INTEGER USEED 

IF (NTOT.EQ@.0) RETURN 

NG = 1 


10 CONTINUE 


12 
15 


.... SUPPRESSION MORTALITY 


IS DEPENDENT UPON SPECIES AND SIZE OF TREE 


IF (PMORT .NE. 9 . AND. KYR .NE. 1) GO TO 15 

WRITE (9,12) 

FORMAT (3X, 3H -1.10X,39H DATA FOR DEAD TREES = SPECIES DBH YEAR) 
CONT INUE 


20 DO 130 Ke=1,NTOT 


50 


60 


70 


80 


I = NCODE(K) 


....CHECK FOR SUPPRESSION 


IF (NOGRO(NG). NE.K) GO TO 130 
NG = NG+1 


....S80RT MORTALITY TYPES 


GO TO (100, 80, 80. BO, 50, 50, 50, 60, 60, 80, BO, 80, 80, 30, BO, 30, 50, 60, 
100, 100,50), I 


....LONG-LIVED EARLY SERAL SPECIES. 


IF (DBH(K). GT. 0. 2#DBHMX(1I)) GO TO 40 

P = 0. 8866+0. 1029/(0. 2-1. /DBHMX (1) )#((DBH(K)-1. )/DBHMX(T) > 
IF (RANF(2).LT.P) 00 TO 130 

GO TO 120 

IF (RANF(2).LT. 0.9895) GO TO 130 

6O TO 120 


.... SHORT-LIVED EARLY SERAL SPECIES. 


IF (RANF(2).LT. 0.628) GO TO 130 
60 TO 120 


....LONG-LIVED MID SERAL SPECIES. 


IF (DBH(K). GT. O. 1#DBHMX(I)) GO TO 70 

P = 0. 96496+0. 01294/(0. 1-1. /DBHMX (1) )#(DBH(K)-1. ) /DBHMX (I) 
IF (RANF(2).LT.P) GO TO 130 

GO TO 120 

P = 0. 97556+0. 02339#DBH(K) /DBHMX (1) 

IF (RANF(2).LT.P) GO TO 130 

60 TO 120 


....LONG-LIVED MID SERAL SPECIES. INCREASING MORTALITY. 


IF (DBH(K). OT. 0. S5*DBHMX(I)) GO TO 90 

P = 0. 8866+0. 1112/(0. 5-1. /DBHMX (1) )#(DBH(K)-1. ) /DBHMX (I) 
IF (RANF(2).LT.P) @O TO 130 

60 TO 120 

P = 1. 0538-0. 1112#DBH(K) /DBHMX (I) 

IF (RANF(2).LT.P) GO TO 130 

@0 TO 120 


.... LATE SERAL SPECIES. 


(0060) 
(0061 ) 
(0062) 
(0063) 
(0064) 
(0065) 
(0066) 
(0067 ) 
(0068) 
(0069) 
(0070) 
(0071) 
(0072) 
(0073) 
(0074) 
(0075) 
(0076) 
(0077) 
(0078) 
(0079) 
(0080) 
(0081) 
(0082) 
(0083) 
(0084) 
(0085) 
(0086) 
(0087 ) 
(0088) 
(0089) 
(0090) 
(0091) 
(0092) 
(0093) 


aaana 


aan 


100 


120 


125 


126 


130 


140 


150 


160 


IF (DBH(K). GT. 0. 1#DBHMX(I)) GO TO 110 

P = 0. 9956+0. 00335/(0. 1-1. /DBHMX (1) )#(DBH(K)—-1. )/DBHMX (1) 

IF (RANF(2).LT.P) @O TO 130 

GO TO 120 

P = 1. 0053-0. 0635#DBH(K) /DBHMX (I) 

IF (RANF(2).LT.P) @O TO 130 

CONT INUE 

NTREES(I) = NTREES(I)-1 
aiicirews CHECK TO SEE IF DEAD TREE CAN STUMP SPROUT. SET KSPRT = -1 
Sulit aie IF TREE CAN SPROUT 

IF (DBH(K). GT. SPRTMN( 1). AND. DBH(K). LT. SPRTMX(I)) KSPRT(I) = -1 

IF (PMORT .NE. 9%) GO TO 126 
FORMAT (3X, 13,F6. 2,14) 

WRITE (9,125) NCODE(K), DBH(K), KYR 


CONTINUE 
DBH(K) = -1.0 
CONTINUE 
|... REWRITE DIAMETERS AND AGES TO ELIMINATE DEAD TREES 
K = 0 


DO 150 I=1,NTOT 
IF (DBH(I).LT.O.) @O TO 150 
K = K+1 
DBH(K) = DBH(I) 
IAGE(K) = IAGE(T) 
NCODE(K) = NCODE(T) 

CONTINUE 

NTOT = K 

CONTINUE 

RETURN 

END 
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20 


30 


70 


SUBROUTINE OUTPUT (KYR, IPLOT) 

COMMON /FOREST/ NTREES(50), DBH( 700), IAGE( 1000), KSPRT (50), NEWTR(50) 
1, SUMLA( 1000), NEW(50), NCODE (700), ISEL (50) 

COMMON /PARAM/ AAA(8, 100), DMAX(50), DMIN( 50), B3(50), B2(50), ITOL(S0) 
1, AGEMX (50), 6(50), SPRTND( S50), SPRTMN( 50), SPRTMX (50), KTIME (50), A1(50) 
2, A2(50), BB( 50), BC(50),C(50), WMIN( 50), IRE@( 50, 4), DBHMX(50), B1(50), 
3B4(50), HMAX (50), DRANGE (50) 

COMMON /PASSX/ A(501, 25), FB(50) 

COMMON /COUNT/ NTOT., NYEAR, IPERT, INT, IN, PMORT 

COMMON /CONST/ NSPEC, SOILG, PMS, TOI, TOII, SD, IGEOG 

DIMENSION BAR(S50) 

REAL LA 

KYR1 = KYR+1 

AREA = 0.0 

TBAR = 0.0 

TOBAR = 0.0 

DO 10 I=1,NSPEC 


FB(I) = O. 
ACKYR1,1) = O. 
BAR(I) = O. 


ACKYR1, 22) = O. 
ACKYR1,23) = O. 
ACKYR1,24) = O. 
ACKYR1, 25)=0. 0 
IF (NTOT.EQ@.0) GO TO SO 
DO 40 J=1,NTOT 
I = NCODE(J) 
IF (I.NE. 6) @0O TO 20 
HT#137. +B2(6)#DBH(J)-B3(6) #DBH( J) ##2 
BAR(I) = BAR(1I)+0. 0242. 09*(HT/100. )*DBH(J)—-0. OO15*(HT/100. )* 
1 DBH(J)##2 
FBIO=3. 2+1. B9#ALOG (DBH( JU) ) 
@0 TO 30 
BAR(I) = BAR(I)+EXP(A1(1)+BB(1I)*#ALOG(DBH(VJ) )) 
LA=(256. *ALOG(DBH( J) )-639. )#2. 5 
FBIO=EXP (A2(1)+BC(1I)#ALOG(DBH( UJ) ) >) 
IF (DBH(J). GE. SO. ) FBIO@LA*(1/C (16) )#DRANGE (I) 
BAR(I) = BAR(I)+FBIO 
ACKYR1, 25) AC KYR1, 25) +FBIO#. 005 
FB(I)=FB(1)+FBIO#. OOS 
GO TO (32, 32, 32, 32, 34, 34, 34, 32, 32, 32, 32, 32, 32, 32, 32, 32, 34, 
1 32, 32, 32, 34), 1 
AREA=AREA + FBIO#C(I)/2.5 
GO TO 40 
AREA=AREA + FBIO#C(I)/2.0 
CONTINUE 
CONTINUE 
DO 60 I=1, NSPEC 
TBAR = TBAR+BAR(I) 
TBAR = TBAR#0O. 005 
DO 70 IV1i=1, NSPEC 
ACAYR1I,1IV1) = ACKYR1I, IVI)+BAR(IV1)#0. OOS 
CONTINUE 
ATOT = NTOT 
Ni = NSPEC+1 
N2 = Ni+1 
N3 = N2+i1 
ACKYRi,N1) = ACKYR1,N1)+TBAR 
TYR(KYR1) = TBAR 
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SUBROUTINE PLOTIN (IPLOT) 
COMMON /FOREST/ NTREES(50), DBH( 700), IAGE (1000), KSPRT(50), NEWTR(50) 
1, SUMLA( 1000), NEW(S0), NCODE (700), ISEL(50) 
COMMON /CONST/ NSPEC, SOILG, PMS, TGI, TGII, SD, IGEOG 
COMMON /COUNT/ NTOT. NYEAR, IPERT. INT, IN, PMORT 
sieanepe INITIALIZE VARIABLES TO START SIMULATION ON BARE PLOT 
Raseeore NTREES CONTAINS NUMBER OF TREES FOR EACH SPECIES 
5) es DBH CONTAINS DIAMETER AT BREAST HEIGHT FOR EACH TREE 
Spee KSPRT IS USED TO FLAG THE TREES THAT CAN SPROUT 
Sauna AREA CONTAINS THE LOCATION OF EACH TREE 
aera loe NOGRO IS USED TO FLAG THE TREES THAT DON"T GROW 
eect IAGE CONTAINS THE AGE FOR EACH TREE 
IPLOT = IPLOT+1 
NTOT = IN 
DO 10 I=1, NSPEC 
IAGE(I) = 0 
NTREES(I) = 0 
10 CONTINUE 
NSPE1 = NSPEC+1 
DO 20 I=NSPE1, 700 
IAGE(I) = 0 
20 CONTINUE 
IF IN NOT EQUAL TO ZERO READ IN TREE SIZES AND SPECIES 
DATA AT END OF DATA FILE SUCDAT 


IF (IN .E@. 0) 6O TO 80 
READ IN THE NUMBER OF TREES OF EACH SPECIES 


READ (5,25) (NTREES(I), I=1, NSPEC) 
25 FORMAT (2113) 
READ (5,30) (NCODE(K), DBH(K), Kei, IN) 
30 FORMAT (7(13, 1X,F6. 2, 7X)) 
40 CONTINUE 
50 CONTINUE 
WRITE (6,60) NTOT 
60 FORMAT (8H NTOT =. I5) 
WRITE (6,70) (NCODE(K), DBH(K), K=1, IN) 


70 FORMAT (1X, 7(15,F6. 2)) 
80 CONTINUE 
RETURN 


END 
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SUBROUTINE PLTREE (KYR, IPLOT) 
COMMON /PARAM/ AAA(S, 100) 
COMMON /FOREST/ NTREES(50), DBH(700), IAGE( 1000), KSPRT(50), NEWTR(50) 
1, SUMLA( 1000), NEW(50), NCODE( 700), ISEL (50) 
COMMON /CONST/ NSPEC, SOILG, PMS. TGI, TGII, SD, IGEOG 
COMMON /COUNT/ NTOT, NYEAR, IPERT, INT, IN. PMORT 
COMMON /PASSX/ A(501, 25), FB(50) 
DIMENSION D(25, 20) 
DATA PI/0. 314159265E01/ 
KYR1 = KYR+1 
BAREA=O. 
DO 5 K=1, NTOT 
BAREA=BAREA+P 1#(DBH(K) /2. )##2 
bo) CONTINUE 
BAREA=(BAREA #* 5. )/10000. 


PRINT OUT TABLE FOR LAST YEAR 
IF (NYEAR. EG. KYR) @O TO 9 


WRITE ONTO TAPES STAND INFORMATION TO BE USED IN APLOT: 
YEAR, STAND BIOMASS. LEAF AREA, STAND BASAL AREA, & 
BIOMASS FOR ABIES AMABILIS, PSEUDOTSUCGA MENZIESII & TSUGA 
HETEOPHYLLA 
RECORDED FOR EVERY OTHER YEAR ONTO TAPES 


IF (MOD(KYR,2) .NE. O) GO TO 100 
WRITE (8,8) KYR, ACKYR1, 22), ACKYR1, 24), BAREA, 
1 ACKYR, 1), ACAYR, 16), ACKYR, 19) 
FORMAT (7(E10. 3, 2X)) 


ONLY WRITE OUT TABLE EVERY SO YEARS 
OR FOR FIRST 5 YEARS 
IF (KYR.LE.5) GO TO 9 

IF (MOD(KYR, 50) .NE. 0) GO TO 100 


WRITE (6,10) KYR, IPLOT. TGI. PMS 
10 FORMAT (1H1, 4HYEAR, 16, 2X, 1L1HPLOT NUMBER, 16, 4X, 6HTGI = .F6. 1, 4X, 
16HPMS = ,F6. 1, 4X, 1BHCEQGRAPHIC REGION: ) 
60 TO (12, 14,16, 18), IGEOG 
12 WRITE (6, 13) 
13 FORMAT (1H+, T86, 22H SOUTH OF SANTIAM PASS) 
60 TO 20 
14 WRITE (6,15) 
15 FORMAT (1H+, T86, 26H SANTIAM PASS TO MT. RAINIER) 
GO TO 20 
16 WRITE (6,17) 
17 FORMAT (1H+, T86, 22H MT. RAINIER TO CANADA) 
60 TO 20 
18 WRITE (6, 19) 
19 FORMAT (1H+, T86é, 18H OLYMPIC PENINSULA) 
BIOM = AC(KYR1, 22) 
20 WRITE (6,22) ACKYR1, 22), A(KYR1, 23), ACKYR1, 24) 
22 FORMAT (1H+, 45H TOTAL FOLIAGE BIOMASS PLUS STEM WOOD BIOMASS 
129H (W/O BARK) TONS PER HECTARE=, F160. 5, 
2/2X, 38H TOTAL NUMBER OF TREES PER 1/5 H PLOT=. Fé. 1, 6X, 4HLEAF 
312H AREA INDEX=, F10. 5) 
WRITE (6,27) BAREA, AC(KYR1I, 25) 
27 FORMAT (1HO, 18H THE BASAL AREA IS. F10.2,14H SQ. M/HECTARE, 6X, 


(0060) 145H TOTAL FOLIAGE BIOMASS (INCLUDES PETIOLES) IS,F8.3//) 


(0061) C 

(0062) C CALCULATE THE NUMBER OF TREES IN EACH DIAMETER CLASS (FOR 
(0063) C EVERY 10 CM). FIRST INITIALIZE D(I.J) TO ZERO. 

(0064) C 

(0065) DO 31 I=1, NSPEC 

(0066) DO 31 Ji, 20 

(0067) D(I, J)20. 0 

(0068) 31 CONTINUE 

(0069) KTOT#A(KYR1, 23) 

(0070) DO 46 K=1, NTOT 

(0071) T=NCODE (K) 

(0072) DC=DBH(K) /10. 

(0073) IF (DC .GT. 20.) WRITE (6,35) DC 

(0074) 35 FORMAT (SH DC= ,F8. 2) 

(0075) IDC=IFIX (DC) 

(0076) IF(IDC .GE. 16) IDC#=16 

(0077) D(I, IDC)=D(1I, IDC) +1. 

(0078) 46 CONTINUE 

(0079) WRITE (6, 47) 

(0080) 47 FORMAT (1H . 6X, 7ZHSPECIES, 13X,. 7H NUMBER, 2X, 6H ABOVE, 1X, 
(0081) +8H FOLIAGE, 

(0082) 120X, 26H NUMBER PER DIAMETER CLASS/24X, 9H OF TREES, 1X, 
(0083) 27H GROUND, 2x, 

(0084) 38H BIOMASS, 1X,40H 10 20 30 40 30 60 70 80 90 100 
(0085) 48H 110 120 | 

(0086) S17H 130 140 150 160+/34X, 8H BIOMASS 

(0087) 611X, 35H-20 -30 -40 -50 -60 -70 -80 -90-100 

(0088 > 725H-110-120-130-140-150-160 ) 

(0089) DO 49 I=1,NSPEC 

(0090) WRITE (6,48) I, (AAACJ, I), J=1,5),NTREES(I), ACKYR, 1), FBI), 
(0091) 1(D(I, J), J=1, 16) 

(0092) 48 FORMAT (1H ,12,1X,5A4 .2X, 16, 2X, F6. 2, 2X, F6. 2, 4X, 16(F 4. 0)? 
(0093) 49 CONTINUE 

(0094) WRITE (6, 50) 

(0095) 50 FORMAT ( ///) 

(0096) C 60 WRITE (6, 65) 

(0097) 65 FORMAT (1HO, 41HSPECIES CODES, DIAMETERS, AND AGES FOLLOW) 
(0098) C WRITE (6,70) (NCGDE(K), DBH(K), IAGE(K), K=1, NTOT) 

(0099) 70 FORMAT (1HO, 7(13, 1X, F6. 2, 1X, 14, 2QH--) /(1HO, 7(13, 1X, Fe. 2, 1X, 14, 2H--) 
(0100) 1)) 

(0101) 100 CONT INVE 

(0102) C 

(0103) C FOR PMORT27, PRINT OUT STAND STRUCTURE AT TIME 700 FOR INPUT TO 
(0104) C NEXT RUN OF CLIMACS 

(0105) C 

(0106) IF (PMORT .NE. 7) GO TO 120 

(0107) IF (KYR .NE. 700) 60 TO 120 

(0108) WRITE (9, 110) (NCODE(K), DBH(K), Kei, NTOT) 

(0109) 110 FORMAT (7(13, 1X, F6. 2, 7X)) 

(0110) 120 CONTINUE 

(O111) C 

(0112) RETURN 


(0113) END 
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FUNCTION RANDOM (NSEED) 
COMMON /SEED/ USEED(13) 
INTEGER USEED 

ISEED = USEED(NSEED) #100. 
DUM = RND(ISEED) 

RANDOM = RND(O) 

RETURN 

END 


SUBROUTINE RANSEE 
COMMON /SEED/ USEED(13) 
INTEGER USEED 


DUM = RND(1) 


DO 10 I = 1,13 
10 USEED(I) = RND(O) 


RETURN 
END 


REAL#4 FUNCTION RANF( IDUM >) 
UNIFORM RANDOM NUMBER GENERATOR 
RANF = RND(O) 


RETURN 
END 


SUBROUTINE SELECT 

COMMON /FOREST/ NTREES(50), DBH( 700), IAGE( 1000), KSPRT (50), NEWTR(50) 
1, SUMLA( 1000), NEW(S0), NCODE (700), ISEL (50) 

COMMON /PARAM/ AAA(B, 100), DMAX(50), DMIN( 50), B3(50), B2(50), ITOL(50) 
1, AGEMX (50), @(50), SPRTND( 50), SPRTMN( 50), SPRTMX(50), ATIME (50), A1(50) 
2, A2(50), BB( 50), BC(50),C( 50), WMIN( 50), IREG( 50, 4), DBHMX(50), B1(50), 
3B4(50), HMAX (50), DRANGE (50) 

COMMON /CONST/ NSPEC, SOILQG, PMS, TCI, TGII, SD, IGEOG 

DO 10 I=1,NSPEC 

ISEL(I) = O 
ore CHECK MAXIMUM WATER STRESS. 
IF (WMIN(I).LT. PMS) GO TO 10 

sted aie CHECK TEMPERATURE REQUIREMENTS. 

IF (TGI.LT. DMIN(I). OR. TOI. GT. DMAX(I)) CO TO 10 

Aeon CHECK SPECIES DISTRIBUTIONS. 

Meee rath REGION 1 = SOUTH OF SANTIAM PASS. 

Sista 6 REGION 2 = SANTIAM PASS TO MT. RAINIER NP. 

Aroereny REGION 3 = MT. RAINIER NP TO CANADIAN BOARDER. 

Ake Yee REGION 4 = OLYMPIC NP. ; 

IF (IREG(I, IGEO@G). NE. 1) @O TO 10 
ISEL(I) = 1 
10 CONTINUE 
RETURN 
END 


(0001) SUBROUTINE SPROUT 


(0002) COMMON /FOREST/ NTREES(S0), DBH( 700), IAGE( 1000), KSPRT(50), NEWTR(50) 
(0003) 1, SUMLA( 1000), NEW(50), NCODE (700), ISEL (50) 

(0004) COMMON /PARAM/ AAA(8, 100), DMAX( 50), DMIN( 50), B3(50), B2(50), ITOL(50) 
(0005) 1, AGEMX (50), 6(50), SPRTND( 50), SPRTMN( 50), SPRTMX (SC), KTIME (50), A1(50) 
(0006) 2, A2(50), BB( 50), BC (50), C( 350), WMIN( 50), IREG(S0, 4), OBHMX(50), B1(50), 
(0007 ) 3B4 (50), HMAX (50), DRANGE (50) 

(0008) COMMON /CONST/ NSPEC, SOILQG. PMS, TGI, TGII, SD, IGEODG 

(0009) COMMON /DEAD/ NOGRO(700) 

(0010) COMMON /COUNT/ NTOT, NYEAR, IPERT, INT. IN, PMORT 

(0011) COMMON /SEED/ USEED(13) 

(0012) INTEGER USEED 

COOKS) Rie Cent tee 

COOV4) ECT ee: SMALLEST AVERAGE STUMP SPROUT iS .1 CM 

COOLS) i Conor: 

(0016) SIZE = .1 

COLOUR ©] aaooe 

COONS) GI ee eee SUM TOTAL NUMBER OF TREES 

(COXON) GS lows 

COOZO) RIG Sea: DETERMINE WHICH SPECIES CAN SPROUT 

COO 21 RC ere ae 

(0022) NW = 0 

(0023) DO 10 I=1, NSPEC 

(0024) IF (SPRTND(I).LE.0O.) GO TO 10 

(0025) IF (KSPRT(I).GE.0) GO TO 10 

(0026) NW = NW+1 

(0027) NEW(NW) = I 

(0028) 10 CONTINUE 

COOLS) CG wea ier ar 

COO3SON ING? ase CHECK FOR SPROUTS 

COO 3)10) sae Cen 

(0032) IF (NW.EQ@.0) GO TO 30 

(OO33) iC Cra 

(OOS4) Cy ae ae CHOOSE RANDOM NUMBER OF SPROUTS 

VOOR) G  seacs 

(0036) NW = NW#RANF(3)+1. 0 

COI) ©] sono 

SOROEISY (GF sa c00 SELECT SPECIES TO SPROUT 

(O03) ia C ate 

(0040) NSPC = NEW(NW) 

COO 4:10) Sar Cree 

COOS20 5G ae SPRTND IS THE TENDENCY FOR THE ITH SPECIES TO STUMP OR 
(0043) Cee ROOT SPROUT. THE VALUE OF SPRTND IS THE AVERAGE NUMBER 
(OO04:4)) Ceo ee OF SPROUTS THAT MIGHT OCCUR WITH A TREE DEATH 

(OO45) EG eee 

SOV) (G  oaacc 

(OOS) ReGen RANDOMLY SELECT NUMBER OF TREES TO SPROUT 

SOOSE)E Ge ase 

(0049) NSPRT = RANF (4) #SPRTND(NSPC) +1 

(0050) DO 20 I=1,NSPRT 

(0051) NTREES(NSPC) = NTREES(NSPC)+1 

(0052) NTOT = NTOT+1 

(0053) IF (NTOT. GT. 700) CALL ERR 

(0054) IAGE(NTOT) = O 

(0055) NCODE(NTOT) = NSPC 

(0056) DBH(NTOT) = SIZE+. 1#(1. O-RANF (5) )##3 

SORE)  (G Ze pac 

(ORORE) (Gg oc STORE DIAMETERS AND AGES FOR NEW SPROUTS 

SOLIS) - Gs los 
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Dale, Virginia H.; Hemstrom, Miles. CLIMACS: A computer model of forest stand development 
for western Oregon and Washington. Res. Pap. PNW-327. Portland, OR: U.S. Department of 
Agriculture, Forest Service, Pacific Northwest Forest and Range Experiment Station; 

1984. 60 p. 


A simulation model for the development of timber stands in the Pacific Northwest is described. 
The model grows individual trees of 21 species in a 0.20-hectare (0.08-acre) forest gap. The 
model provides a means of assimilating existing information, indicates where knowledge is 
deficient, suggests where the forest system is most sensitive, and provides a first testing ground 
for hypotheses. Model verification simulations are included for up to 500 years on various sites. 


Fire, wind, or clearcutting can occur at intervals and intensities specified by users. The model 
was developed by modifying an existing forest succession simulator of eastern deciduous for~- 
ests. Birth, growth, and death of individual trees are functions of existing light and temperature 
conditions, competition and species characteristics. Modifications of the existing simulator in- 
clude tree height growth being related to temperature and moisture conditions, the foliage bio- 
mass to diameter relationship being more realistic, and five mortality classes and shade tolerance 
classes being defined. 


Keywords: Succession, models, simulation, community dynamics (plant). 
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